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ABSTRACT. Squamate s ee urds, snakes, and their kin 
such as amphisbaenians, or “worm lizards”) represent 
the world’s most diverse Ade of terrestrial vertebrates 
with ~11,000 described extant species, representing 
key components in many of the world’s most diverse 
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ecosystems. With an evolutionary history dating back 
at least to the Middle Triassic at 242 Ma, the squamate 
Tree of Life also features numerous diverse but extinct 
branches, with hundreds of fossil species found all 
over the world. Despite their biological relevance both 
today and in the geological past, there remains a 
centuries-old controversy on how the major lineages of 
squamates are related to each other, with a Guee 
impact on studies in ecology, evolution, paleontology, 
toxinology, and other fields. Here, we provide a a 
itor. overview of this long research tradition, 
from 19th century naturalists to 21st century phylo- 
genomics, with special emphasis on sever ral recent 
advances over the last two decades. These insights 
have had a dramatic effect on our understanding of the 
squamate Tree of Life and clarify several possible 
future research agendas. We provide an integrative 
perspective deriv ad from genomics, morphology, and 
the fossil record and propose several points of 
synthesis in our current knowledge of broadscale 
squamate evolution and systematics. Key topics of 
interest include dating the origin and early evolution 
of lizards, the phylogenetic origin of snakes, the 
evolution of venom, recent agreements between 
morphological and molecular squamate evolutionary 
trees, genomic patterns of evolution, and the integra- 
tion th morphological and molecular data sets. We 
conclude by providing perspectives on possible 
advancements in the field, directing researchers to 
promising future lines of investigation that are 
necessary to further expand our synthetic knowledge 
of squamate evolution, 


Key words: Le pidosi wiria, Squamata, Lizards, 
Snakes, Amphisbaenians, Phylogeny, Evolution, Phy- 
logenomics, Morphology, C Yombined evidence, Sys- 
tematics, Total evidence dating, Fossils 


INTRODUCTION 


Squamates (lizards, snakes and amphis- 
baenians) are the most diverse known clade 
of extant terrestrial vertebrates with c: 
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11,000 species, nearly twice as that of 
mammals (Uetz and HoSek, 2021). With 
an evolutionary history spanning more than 
242 million years (Simdes et al., 2018), 
squamates have inhabited nearly all terres- 
trial and aquatic environments from the 
distant past to today, including all of the 
world’s oceans during the Cretaceous, and 
from Antarctica to the Arctic Circle during 
warmer periods (Pianka and Vitt, 2003; Vitt 
and Caldwell, 2013). An accurate under- 
standing of the squamate Tree of Life is 
thus a crucial basis for reconstructing the 
evolutionary history of a considerable com- 
ponent of the world’s biodiversity. A robust 
phylogenetic tree is therefore the necessary 
starting point towards a broader under- 
standing of the evolvability of functional 
traits (e.g., Wiens et al., 2006), distribution 
of extinction risk (e.g., Tonini et al., 2016), 
origin of life history strategies (e.g., Pyron 
and Burbrink, 2014), and myriad other 
macroevolutionary dynamics (e.g., Title 
and Rabosky, 2017). 

Although we are now in the fourth 
decade of quantitative phylogenetic analyses 
of squamates (beginning with Estes et al., 
1988; see below) and their various sub- 
groups, numerous questions remain unan- 
swered. Morphological data sets differ 
strongly in the placement of key fossil 
groups such as Mosasauria and Paramacel- 
lodidae (see Pyron, 2017) and relationships 
among the major extant families. Molecular 
data sets differ in the placement of key 
groups such as Dibamidae; the relationships 
of Iguania, Anguimorpha, and Serpentes; 
and relationships within Iguania (see Bur- 
brink et al., 2020). Perhaps most famously, 
many morphological and molecular data 
sets appear fundamentally irreconcilable 
with respect to the placement of taxa such 
as Iguania (Losos et al., 2012). Although 
recent genomic and morphological data sets 
resolve some of these issues (Sim6es et al., 
2018; Singhal et al., 2021), numerous 
branches are still unresolved. These linger- 
ing questions about the squamate Tree of 


Life have affected our ability to answer 
fundamental questions about squamate 
systematics and evolutionary biology. 

First, what is the best phylogenetic 
hypothesis for Squamata, including both 
extinct and extant lineages? Are the various 
traditional morphological hypotheses driven 
by widespread homoplasy (e.g., Losos et al., 
2012)? Are evolutionary processes such as 
paralogy, incomplete lineage sorting, and 
horizontal gene transfer or hybridization 
driving conflicting topologies in genome- 
scale data (e.g., Siu-Ting et al., 2019)? 
Could other sources of morphological 
(e.g., soft tissue; Siegel et al., 2011) or 
molecular (e.g., synteny mapping; Sevillya 
and Snir, 2019) data offer alternative 
perspectives on species relationships? 
Could both morphological and molecular 
estimates be confounded by analytical 
problems such as gene tree error (Roch 
and Warnow, 2015) or character construc- 
tion and coding (SimGes et al., 2017b)? Can 
refinement of characters (Simdes et al., 
2018) and integration of ontogeny (Wiens 
et al., 2005) help to alleviate the apparent 
molecular-morphological discordance? 

Second, once a higher level resolution of 
the squamate Tree of Life is achieved, what 
can we learn about the evolutionary history 
of snakes, lizards, and worm-lizards? What 
are the earliest evolving lineages of squa- 
mates and their closest reptile relatives? 
When did crown squamates originate? Do 
snakes have a marine, terrestrial, or fossorial 
ancestry (Vidal and Hedges, 2004; Lee, 
2005a)? Did venom originate only once or 
multiple times across lizards and_ snakes 
(e.g., Fry et al., 2006)? The questions that 
such an estimate can help answer are nearly 
endless. 

Here, we review the recent history and 
future directions of research into the 
squamate Tree of Life. We detail many of 
the most contentious relationships, the data 
and methods used to arrive at various 
results, and the conclusions reached by 
recent authors regarding phylogenetic to- 
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pologies and evolutionary histories in the 
group. We then offer a synthesis and 
perspective on where knowledge currently 
stands for these contentious branches. 
Finally, we detail potentially fruitful ave- 
nues for future investigations. We believe 
that promise exists in both data and 
methods to drastically increase our knowl- 
edge of lizard and snake relationships and 
that these will provide grounds for an 
increasingly broadscale synthetic view of 
squamate evolution. 


A BRIEF OVERVIEW OF SQUAMATE 
DIVERSITY 


Extant Lizards and Snakes 


The phylogenetic diversity and taxonomic 
history of extant squamates (as well as 
extinct lineages; see below) is staggeringly 
complex, with a bewildering array of clade 
names, definitions, and transfers of lineages 
between taxa as phylogenetic understanding 
increased (see Conrad, 2008; Vitt and 
Caldwell, 2013). Recounting these names 
and their evolving composition over time is 
too substantial an undertaking to present 
here. Instead, we will briefly outline the 
major clade names and their composition 
from the recent study of Burbrink et al. 
(2020), who analyzed higher level relation- 
ships within Squamata with genome-scale 
data. Those authors provide a_ baseline 
understanding of phylogeny and taxonomy 
for a series of higher level clades for which 
consensus has emerged since the formaliza- 
tion of the “molecular” hypothesis of 
squamate relationships (see Vidal and 
Hedges, 2005, 2009; Pyron et al. 2013). 
Many of these taxa (e.g., Serpentes, the 
snakes) have been recognized in some form 
or another since initial historical classifica- 
tions were first erected. We also comment 
on the diversity and composition of each 
group (see Vitt and Caldwell, 2013; Uetz 
and HoSek, 2021). 

Although species, genera, and even entire 
families and subfamilies have been moved 


into or out of these taxa based on morpho- 
logical and molecular analyses over the 
years, we present them here as defined by 
most recent molecular phylogenies (e.g., 
Vidal and Hedges, 2005; Pyron et al., 2013; 
Burbrink et al., 2020). Serpentes is perhaps 
the most well characterized, as few extant 
snakes have ever been misclassified as 
lizards or vice versa. This group was 
historically classified (see below) into the 
blindsnakes and threadsnakes (Scolecophi- 
dia); boas, pythons, and allies like pipe- 
snakes, sunbeam snakes, and _ shieldtails 
(Henophidia); and advanced snakes (Cae- 
nophidia), comprising many of the most 
common species worldwide and all medi- 
cally significant venomous species. Scoleco- 
phidia and Henophidia are most likely not 
monophyletic, although Caenophidia and 
Caenophidia+“Henophidia” (Alethinophi- 
dia) are strongly supported in most recent 
phylogenies (see Burbrink et al., 2020; 
Singhal et al., 2021). Snakes represent 
~3,900 species (Uetz and HoSek, 2021) 
occurring on every continent except Ant- 
arctica and in nearly every terrestrial and 
aquatic environment, spanning a massive 
range of morphological and ecological 
diversity. 

Iguania is a diverse clade of ~2,000 
species of primarily diurnal, terrestrial and 
arboreal lizards including charismatic and 
well-known taxa such as Anolis, chameleons 
and agamas, iguanas, and other diverse 
assemblages. Whereas most iguanian fami- 
lies (Pleurodonta) inhabit Neotropical re- 
gions in the Americas, acrodontans (agamids 
and chameleons) inhabit Old World conti- 
nents. The exception to that are some 
pleurodont families with species in Mada- 
gascar (Opluridae) and the Galapagos and 
Fiji (iguanidae), which give puzzling indi- 
cators of either Gondwanan vicariant origin 
or long-distance oceanic dispersal (Noonan 
and Chippindale, 2006; Noonan and Sites, 
2010). Surprisingly, the most abundant 
fossil record of extinct pleurodontan igua- 
nians occurs in East Asia (e.g., Gao and 
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Norell, 2000) whereas the oldest known 
acrodontans came from what is today South 
America and Africa (Sim6ées et al., 2015b; 
Apesteguia et al., 2016), further confound- 
ing our understanding of their biogeograph- 
ic history. 

Several apparently “primitive” or plesio- 
morphic morphological traits, such as a 
fleshy, lobular tongue used for lingual prey 
prehension, led historical authors and early 
morphological systematists to place Iguania 
as the earliest diverging lineage of squamate 
(see below; Losos et al., 2012). Relation- 
ships among the many pleurodontan igua- 
nian families have proven intractable to 
estimate even from genome-scale data sets, 
apparently because of a rapid evolutionary 
radiation with short intervening branch 
lengths (Townsend et al., 2011; Singhal et 
ales, 2020): 

Another group of widely recognized 
lizards has often been referred to as 
Anguimorpha (see below), but was referred 
to as Anguiformes (Conrad, 2006) by 
Burbrink et al. (2020). This includes the 
famous monitors (Varanidae) such as the 
Komodo dragon; the venomous Gila mon- 
ster and beaded lizards (Helodermatidae); 
the Mexican knob-scaled lizards (Xenosaur- 
idae); the Chinese crocodile lizard (Shini- 
sauridae); the American legless lizards 
(Anniellidae); the galliwasps (Diploglossi- 
dae); and the alligator lizards, glass lizards, 
and slow worms (Anguidae). Despite a 
relatively small number of species (~250) 
compared with other diverse squamate 
clades, anguiforms nonetheless contain a 
number of exceptionally distinct species 
that are ingrained in the public conscious- 
ness. The group is primarily tropical in 
distribution, with only a few species in 
temperate regions, but occurring on every 
continent except Antarctica. Some of the 
oldest known records of the group come 
from the Late Cretaceous of North America 
and the Early Cretaceous of the UK— 
represented by Dorsetisaurus (Hoffstetter, 
1967; Prothero and Estes, 1980). Anguids in 


particular become more abundant during 
the Eocene in North America, including the 
now extinct and heavily armored Glypto- 
saurinae (Sullivan, 1979; Estes, 1983). 
Scincoidea or Scinciformata (see Vidal 
and Hedges 2005, 2009; Pyron et al. 2013) 
represents a distinct clade of terrestrial, 
arboreal, or fossorial lizards, typically pre- 
senting hard shiny dorsal scales or forms of 
plate armor, in the families Scincidae 
(skinks), Gerrhosauridae (plated lizards), 
Cordylidae (girdled lizards), and Xantusii- 
dae (night lizards). Gerrhosaurids, cordyl- 
ids, and xantusiids contain only a few dozen 
species each and are restricted to southern 
Africa (Cordylidae), Africa and Madagascar 
(Gerrhosauridae), or the desert southwest 
of North America, Cuba, and _ tropical 
Central America (Xantusiidae). In contrast, 
skinks are hyperdiverse (~1,700 species) 
and distributed worldwide in nearly every 
region and community containing squa- 
mates, with exceptional ecomorphological 
diversity and habitat variability. Fossils of 
scincoids are relatively rare in the fossil 
record, but early scincoids (or close rela- 
tives) may comprise one of the first 
squamate clades to show up in the fossil 
record (see below—Paramacellodidae). 
Lacertoids or Laterata represents a di- 
verse group (~1,000 species) of primarily 
terrestrial, diurnal, highly active, and con- 
spicuous lizards, along with an enigmatic 
clade of burrowing, limb-reduced species. 
Lacertidae (wall or true lizards) are a 
common feature of squamate communities 
across Europe, Africa, the Middle East, and 
Asia. Gymnophthalmidae (spectacled _liz- 
ards) are ecologically widespread in Central 
and South America, including savannas, 
deserts, and montane and lowland rain- 
forests, and are typically smaller, more 
nocturnal, and more cryptic or fossorial 
than the other terrestrial lacertoids. Teiids 
are a widespread component of typically 
arid- or semiarid-adapted squamate assem- 
blages across North and South America, but 
including some rainforest and semiaquatic 
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species, ranging from small (~10 em total 
length) to large (>1 m) body sizes. 

Nested within Lacertoidea as the likely 
sister lineage to Lacertidae (see below), are 
the enigmatic Amphisbaenia, the “worm 
lizards.” These are limb-reduced or limbless 
species, usually with stout bodies and heads 
adapted for a fossorial burrowing lifestyle, 
with a highly fragmented and disjunct 
distribution across Florida, Central Amer- 
ica, the Caribbean, northern South Amer- 
ica, sub-Saharan Africa, the Iberian 
Peninsula and adjacent northern Africa, 
Madagascar, and the Arabian Peninsula 
(Vitt and Caldwell, 2013). Their phyloge- 
netic placement has been contentious across 
various morphological and molecular analy- 
ses through time (see below). The early 
fossil record of amphisbaenians includes 
extremely well preserved Eocene and Oli- 
gocene rhineurids within the genus Spatho- 
rhynchus (Berman, 1973, 1977; Miiller et 
al., 2016). However, the fossil record of 
amphisbaenians before the Eocene includes 
several purported intermediate forms that 
have subsequently been reassigned to other 
clades by morphological phylogenetic anal- 
yses, such as Sineoamphisbaena (Wu et al., 
1993; Kearney, 2003; Gauthier et al., 2012) 
and Slavoia (Conrad, 2008; Tatanda, 2016; 
Simoes et al., 2018), both from the Late 
Cretaceous of Mongolia. 

Gekkota represents a diverse lineage of 
highly distinct and well-known species, 
consisting of more than 2,000 species. 
Gekkotans are classified into several extant 
families (Diplodactylidae, Carphodactyli- 
dae, Eublepharidae, Sphaerodactylidae, 
Phyllodactylidae, and Gekkonidae), typically 
distinguished from all other squamates by 
their reduced irregular granular scales and 
expanded toe pads that allow many species 
to scale smooth vertical surfaces and even 
hang upside down by Van Der Waals forces 
generated by subdigital lamellae (Autumn et 
al., 2002). Females typically only lay 1-2 
eggs per clutch, in contrast to the wide 
range of fecundity observed in other squa- 


mate lineages. Gecko species are found 
worldwide, with their highest diversities in 
tropical and warm desert regions. Nested 
within geckos are the flap-footed lizards 
(Pygopodidae), a small clade (~35 species) 
of limb-reduced or limbless forms of highly 
specialized geckos found in Australia and 
New Guinea. 

In both morphological and molecular 
phylogenies, Gekkota is typically estimated 
as an early-diverging lineage, subject only to 
the placement of Iguania and Dibamia. The 
gekkotan fossil record is relatively scarce, 
with the most informative fossils coming 
from the Late Cretaceous of Mongolia 
(Daza et al., 2014), but it includes the 
oldest adhesive toepads preserved in amber, 
indicating this remarkable adaptation had 
already evolved at least by 100 Ma (Arnold 
and Poinar, 2008). Stem gekkotans also 
display adaptations for scansoriality, despite 
not yet having adhesive toepads on the basis 
of their hand morphology (Daza et al., 2014; 
Simoes et al., 2017a). 

Finally, Dibamia (Dibamidae; blind 
skinks) is a small group (~25 species) of 
highly derived limb-reduced or limbless 
fossorial species found primarily in mesic 
temperate or tropical habitats from South- 
east Asia to western New Guinea (Diba- 
mus), with a single species (Anelytropsis) in 
the deserts of northeastern Mexico. They 
are primarily insectivorous, with reduced 
eyes, and typically lay a single egg with a 
hard calcified shell, unlike the leathery eggs 
of most other squamates. Their reduced 
body form, with convergence on a burrow- 
ing body shape and associated adaptations 
in skull morphology, makes their phyloge- 
netic affiliations uncertain by morphological 
data (see below). Most molecular analyses 
place them as an early-branching lineage of 
Squamata, either the earliest diverging 
lineage, the sister lineage of Gekkota (with 
Gekkota+Dibamia forming the earliest di- 
verging squamate lineage), or as the sister 
lineage to all other squamates excluding 
Gekkota. As with Iguania, rapid evolution- 
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ary radiations with short internodes and 
long descendant branches make decisive 
resolution of these early squamate nodes 
difficult, even with genome-scale data (Sin- 
ghal et al., 2021). The complete absence of a 
fossil record of dibamids makes the matter 
of understanding the origin of dibamids 
even more complex, because it prevents 
access to the early morphology of the group 
and calibration points for its origin. 


Key Fossil Taxa 


Fossils are known throughout the entire 
evolutionary history of Squamata, from early 
stem taxa in the Triassic (Figs. 1A—C) to 
Pleistocene ancestors of extant species. (For 
broad previous reviews, see Estes, 1983; 
Evans, 2003; Caldwell, 2019.) Many extinct 
species are easily assignable to crown 
lineages such as Serpentes (snakes) (Figs. 
1D, E) or Gekkota (geckos) and provide 
clear insight into the evolutionary history of 
those groups (Daza et al., 2014; Caldwell, 
2019). However, several diverse enigmatic 
groups of extinct species form entirely 
distinct evolutionary branches within crown 
squamates, but of uncertain phylogenetic 
affinity to extant lineages. Among the 
several fossil groups of squamates, we here 
provide a stronger focus on three clades that 
have had a large effect on reconstructing 
the squamate Tree of Life, either because of 
their age or taxonomic richness: Paramacel- 
lodidae, Borioteiioidea (or Polyglyphano- 
dontidae), and Mosasauria. 

Paramacellodids (Figs. 1F—I) are perhaps 
the least understood groups of fossil squa- 
mates. They primarily appear to be terres- 
trial lizards but have an extremely 
incomplete fossil record, mostly composed 
of fragmentary and disarticulated fossil 
remains that complicates any further infer- 
ence of ecomorphology. However, this poor 
fossil record is a consequence of the age of 
the group: They were mostly abundant 
between the Middle Jurassic and the end 
of the Early Cretaceous (Bittencourt et al., 


9020), an age from which the terrestrial 
vertebrate fossil record is scarce globally 
(Close et al., 2019) and furthered biased 
against small-bodied vertebrates (Evans, 
2003). On the other hand, the age of 
paramacellodids makes them an incredibly 
important group, because they represent 
the oldest recognizable clade in the evolu- 
tionary history of squamates to diversify and 
occupy multiple continents (Bittencourt et 
al., 2020). Their fossils have been found in 
the Middle Jurassic of Europe and possibly 
Morocco; the Late Jurassic of North Amer- 
ica, Europe, Africa, Central and East Asia 
(China); the Early Cretaceous of Europe, 
Africa, North America, South America, 
Central Asia (Siberia), and East Asia (Japan 
and Mongolia); and the Late Cretaceous of 
Europe (Waldman and Evans, 1994; N ydam 
and Cifelli, 2002; Folie and Codrea, 2005; 
Haddoumi et al., 2016; Bittencourt et al., 
2020). 

Owing to the fragmentary nature of most 
fossils, the primary diagnostic feature unit- 
ing paramacellodids lies in their dentition: 
teeth labiolingually expanded at their bases, 
along with lingually concave tooth apices 
(Figs. 1G-—I). These characteristics for 
dental morphologies are relatively rare 
within squamates and only found in con- 
junction in paramacellodids. Another highly 
relevant phylogenetic character from this 
lineage is the frequent presence of body 
osteoderms (Fig. IF). This osteological 
feature has been one of the most important 
characters for placing paramacellodids as 
closely allied to scincoids in most phyloge- 
nies (see below). It is likely that the 
discovery of more complete specimens of 
paramacellodids, or the redescription of 
poorly known but fairly complete taxa 
(e.g., Sharovisaurus; Hecht and Hecht, 
1984), will be key to illuminate the place- 
ment of this group in the squamate Tree of 
Life. 

From the Late Cretaceous, another 
taxonomically and geographically diverse 
group of squamates has drawn a lot of 
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Figure 1. Sample of extinct squamate diversity. (A-C) Holotype of the oldest stem squamate Megachirella wachtleri (PZO 628), 
from the Middle Triassic of Italy. (A) Whole skeleton in dorsal view. (B) Reconstructed skull in dorsal view. (C) Computed 
tomography—scanned and segmented skull in dorsal view. (D, E) The early legged snake Najash rionegrina, from the Late 
Cretaceous of Rio Negro Province, Argentina. (D) Computed tomography—scanned and segmented skull in right lateral view 
(MPCA-PV-500). (E) Hind limbs and cloacal region (MPCA-PV-400). (F—H) Polyglyphanodon sternbergi (Borioteiioidea), from the 
Late Cretaceous of Utah, USA. (F) Postcranium and osteoderms of Paramacellodus oweni (Paramacellodidae—NHMUK R8209) 
from the Early Cretaceous of the UK. (G—l) Holotype of Neokotus sanfranciscanus (Paramacellodidae—IGC-P 0085), from the 
Early Cretaceous of Minas Gerais, Brazil. (G) Maxillary dentition in occlusal view. (H) Detail of maxillary tooth. (I) Dentary in 
occlusal view. (J) CM 9188 in dorsal view. (K) Holotype skull of P. sternbergi (NMNH 15477). (L) Reconstructed skull anatomy of P. 
sternbergi. (M, N) Skull of Tylosaurus nepaeolicus (Mosasauria—former holotype of 7. ‘kansasensis—FHSM VP 2295). (M) 
Articulated skull in left lateral view. (N) Skull anatomy in left lateral view. Abbreviations: Art-Pra, articular+prearticular; Ang, angular; 
At, Atlas; Ax, Axis; Boc, basioccipital; Bsp, basipterygoid; C, coronoid; Ca.V., caudal vertebrae; CB, compound bone; D, dentary; 
D.V., dorsal vertebrae; Ect, ectopterygoid; F, frontal; Fe, femur; il, lium; J, jugal, L, lacrimal; La.W., labial wall of the dentary; Li.Co., 
lingual concavity; M, maxilla; N, nasal; Ost, osteoderms; P, parietal; Pal, palatine; PM, premaxilla; Po, postorbital; PoF, postfrontal; 
POF, postorbitofrontal; PrF; prefrontal; Pro, prootic; P.St.Pr., parietal supratemporal processes; Ptg, pterygoid; Q, quadrate; RAP, 
retroarticular process; San, surangular; Spm, septomaxilla; Soc, supraoccipital; Spl, splenial; Sq, squamosal; St, supratemporal; 
S.V., sacral vertebrae; Ti, tibia; T.Ba., tooth base; Vo, vomer. Scale bars = 50 jm (H), 500 um (G, I), 5 mm (B-D), 10 mm (A, E, F, 


K, L), and 100 mm (J, M, N). 
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attention for reconstructing the squamate 
Tree of Life. This group has been named 
and rediagnosed several times in recent 
decades, including Polyglyphanodontinae 
(Estes, 1983), Borioteiioidea (Nydam et 
al., 2007), Polyglyphanodontidae (Conrad, 
2008), and Polyglyphanodontia (Gauthier et 
al., 2012) (Figs. 1J—L), as more species have 
been described and anatomical knowledge 
has increased. Historically, borioteiioids 
have typically been linked to teiioids, 
initially as a subfamily of Teiidae (Estes, 
1983) on the basis of gross morphological 
similarity. Later estimates from quantitative 
phylogenetic analyses place the group as the 
sister lineage to Teiidae (Conrad, 2008) or 
Teiioidea (Nydam et al., 2007; Simoes et al., 
2018) on the basis of several coded similar- 
ities between their dentition and that of 
teiioids, such as the heavy deposits of 
cementum at tooth bases, the presence of 
labiolingually expanded teeth, and the 
semicircular shape of the resorption pits 
(Nydam et al., 2007; but see below for a 
distinct hypothesis proposed by Gauthier et 
al., 2012). 

Borioteiioids were geographically wide- 
spread throughout North America and East 
Asia during the Late Cretaceous and may be 
an ancestral lineage involved in the dispers- 
al of teiioids to South America (Nydam et 
al., 2007). Their known fossils include 
individuals with a wide variety of dental 
morphologies, including species with a 
molarlike multicuspid dentition resembling 
that of mammals, as found in Peneteius 
(Nydam et al., 2000). Another example of 
extreme adaptation in the group is found in 
Polyglyphanodon sternbergi, by far the most 
completely known borioteiioid from North 
America with nearly perfectly articulated 
individuals (Figs. 1], K), including an 
ontogenetic series and possible evidence of 
sexual dimorphism (Simdes et al., 2016). 
Polyglyphanodon had an extreme level of 
labiolingual expansion of its dentition, with 
tiny denticles on the occlusal surface that 
were apparently highly adapted for crop- 


ping vegetation (Nydam and Cifelli, 2005). 
Polyglyphanodon and the Chinese borio- 
teiioid Tianyusaurus further possessed an 
extremely unusual feature among squa- 
mates—the presence of a complete lower 
temporal bar forming a double temporal 
fenestration of the skull (Figs. 1], 2, char. 3). 
This complete lower temporal bar is a 
condition reminiscent of early diapsid rep- 
tiles, such as archosaurs, that was reac- 
quired at least twice independently within 
borioteiioids, but never in any known other 
group of living or extinct squamate (Simoes 
et al., 2016). Owing to the completeness of 
several specimens of borioteiioids from 
North America and Asia, particularly Poly- 
glyphanodon and Gilmoreteius, borio- 
teiioids represent one of the very few 
clades of extinct squamates for which their 
morphology is extremely well characterized 
and generally well understood. 

The third and final group of interest for 
the discussions below is Mosasauria (Figs. 
1M, N), a diverse group (~80 valid species) 
of mostly marine lizards with an extraordi- 
nary range in size, from ~5 cm in length in 
the dolichosaurid Aphanizocnemus libanen- 
sis (Dal Sasso and Pinna, 1997) up to an 
estimated 17 m in Mosasaurus hoffmannii 
(Lingham-Soliar, 1995)—a size range span- 
ning at least two orders of magnitude. 
Mosasaurians include early evolving lineag- 
es with only a few described species, such as 
Dolichosauridae (e.g., Dolichosaurus and 
Pontosaurus) and Aigialosauridae (e.g., Ai- 
gialosaurus) (DeBraga and Carroll, 1993: 
Campbell Mekarski et al., 2019). These taxa 
are generally considered to be the forerun- 
ners of mosasaurs—the giant aquatic lizards 
that occupied all of the world’s oceans by 
the end of the Cretaceous (Poleyn et al., 
2014; Jiménez-Huidobro et al., 2017). 
Mosasaurs (Mosasauridae) have traditional- 
ly been subdivided into two main groups, 
the Mosasaurinae and Russellosaurina (Bell, 
1997). The first phylogenetic analysis of the 
group using Bayesian inference ‘suggested 
the latter might not be monophyletic 
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Figure 2. Sample of morphological characters historically used for phylogenetic analyses of squamates. (A) Skull drawing of 
Sphenodon punctatus (Sphenodontia—from FMNH 11118) in left lateral view. (B) Skull drawing of Polychrus marmoratus (\guania, 
Polychrotidae—AMNH R 148543) in left lateral view. (C) Skull of S. punctatus in dorsal view (MCZ R4702). (D) Skull of Varanus 
salvator (Anguimorpha, Varanidae—TMP 1990.007.0031) in dorsal view. (E) Skull of S. punctatus in anterior view (FMNH 11113) 
(F) Skull of Iguana iguana (Iguania, Iguanidae—UAMZ, uncatalogued specimen) in anterior view. (G) Braincase of V. salvator in 
left lateral view (TMP 1990.007.0031). (H) Braincase of Broadleysaurus major (Scincoidea, Gerrhosauridae—AMNH 173621) in 
right ventrolateral view. (I) Premaxilla and anterior palate of Anilius scytale (Serpentes, Anilioidea—MCZ R-19537) in ventral view 
Pectoral girdle of B. major in ventral view (AMNH 173621). Abbreviations: Cl, clavicle; Co, coracoid; Cr.Int., crista interfenestralis; 
Cr.Tb., crista tuberalis; F, frontals; Fn.Ov., fenestra ovalis; Icl, interclavicle; J, jugal; LARST, lateral aperture of recessus scalae 
tympani; P, parietal; Po, postorbital; PoF, postfrontal; Pst, presternum; Sca, scapula; Sq, squamosal; X—XIl, cranial nerves X—XIl 
Phylogenetic character abbreviations: 1, squamosal, dorsal process; 2, quadratojugal; 3, jugal, posteroventral process; 4 
postfrontal; 5, prootic, alar crest; 6, parietals, fusion to each other; 7, premaxillae, fusion to each other; 8, crista tuberalis (division of 
the metotic fenestra); 9, premaxilla, medial dentition; 10, coracoid, anterior emargination; 11, clavicle, secondary anterior 
curvature. Character states within parentheses indicate absence (0), or presence (1) 


(Simoes et al., 2017d), but it was later least two independent land-to-water evolu- 
sustained by relaxed clocks (Madzia and tionary transitions (Caldwell and Palci, 
Cau, 2020). Interestingly, several recent 2007; Simoes et al., 2017d). 

phylogenetic analyses of mosasaur relation- Mosasaurians have played a central role 
ships suggest an independent acquisition of in controversies over morphological analysis 
aquatic adaptations in the group, with at of squamate relationships, specifically owing 
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to their proposed sister group relationship 
to snakes, or as an early branching lineage of 
squamates (see below). Of special interest 
to broader squamate relationships is a 
better understanding of the morphology 
and phylogeny of early mosasaurians that, 
historically, have been far less studied than 
mosasaurids. Some key taxa in that regard 
are the dolichosaurids Pontosaurus and 
Adriosaurus (Lee and Caldwell, 2000; 
Pierce and Caldwell, 2004; Caldwell, 
2006), which along with the aigialosaurid 
Aigialosaurus (DeBraga and Carroll, 1993; 
Caldwell and Dutchak, 2006) represent the 
best preserved early mosasaurians. Deeper 
insights into those specimens from high- 
resolution computed tomography scanning, 
hopefully supplemented by discoveries of 
new well-preserved specimens, will play a 
key role in resolving the placement of 
mosasaurs. 


MORPHOLOGICAL HYPOTHESES 
The Precladistic Era 


Attempts to understand the squamate 
Tree of Life have a long and complex 
history, including simple pre-evolutionary 
classifications in the 18th century from gross 
anatomy, complex morphologically based 
taxonomies informed by evolutionary theory 
in the 19th century, proto-cladistic graphical 
representations of phylogeny in the early 
20th century, and statistical inferences in 
recent decades from morphological and 
molecular data sets. In the pre-evolutionary 
era, we can mark the start of squamate 
systematics (as with most animal groups) by 
the work of Linnaeus (1758), with Laurenti 
(1768), Daudin-(1802), Oppel CISD) Gtv= 
ier (1817), Duméril and Bibron (1839), 
Fitzinger (1843), and Gray (1845) expand- 
ing our knowledge of species richness and 
erecting increasingly sophisticated classifi- 
cations on the basis primarily of morpho- 
logical features, both internal and external. 
In the evolutionary era, Cope (1864), 
Boulenger (1885), Nopesa (1923), and many 


others began to expand and refine these 
classifications in an evolutionary and quasi- 
phylogenetic approach. 

These early works are too numerous and 
expansive in scope, and_ their varying 
taxonomic concepts and nomenclatural 
changes too complex to summarize easily 
here. Yet, from them, we derive many of the 
clade names and natural groupings that are 
still recognized today. The extremely de- 
tailed anatomical investigations of these 
early authors quickly led to the recognition 
of robustly supported “clades” defined by 
morphological “synapomorphies,” even if 
these concepts were not yet developed or 
fixed in the minds of investigators at the 
time. Thus, many genera, subfamilies, and 
families were erected on the basis of clear 
anatomical affinities in the 18th and 19th 
centuries, and many of these hypotheses 
have stood the test of time as both 
morphological and molecular systematics 
have increased in precision and accuracy. 
In contrast, processes such as ecomorpho- 
logical convergence, particularly for limb- 
reduced and burrowing body forms, made 
the early allocation of many such taxa 
difficult and uncertain. Similarly, hyper- 
diverse groups without clearly hierarchical 
diagnostic morphological characters (2 
colubroid snakes, geckos, and_ skinks) 
proved challenging to classify in a coherent 
way below the “family” or “superfamily” 
level. Thus, compelling resolution of many 
branches has only been possible in recent 
decades with molecular data, yielding many 
phylogenetic novelties (see Pyron et al., 
2013). Similarly, clade names such as 
“Sauria” for lizards persist even in the 
present day, despite the fact that Sauria has 
long been known to be paraphyletic with 
respect to Serpentes and Amphisbaenia:; 
both snakes and worm lizards are simply 
derived lizards. 

Two major syntheses in the 20th century 
(Camp, 1923; Underwood, 1967) represent 
the first systematic attempts to unify avail- 
able anatomical data and classify lizards and 
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snakes in a comparative framework with 
explicit notions of common ancestry and 
ancestor-descendant relationships among 
groups. Whereas the latter postdates the 
English edition of Hennig (1966) by one 
year, both are precladistic in method and 
represent the culmination of premodern 
thought in squamate systematics on the 
basis of morphological examination. Taxa 
were grouped and arranged by an implicitly 
subjective understanding of shared charac- 
ters, but generally without formal recogni- 
tion of this principle or quantitative analysis 
thereof. Camp (1923) is by far the most 


influential of all of the precladistic phyloge- 


netic studies, depicting several aspects of 
squamate relationships that would dominate 
morphology-based research in squamate 
evolution for the next century. For instance, 
he proposed an early dichotomy between a 
group formed by iguanians and geckos 
(Ascalabota) and all other squamates, fol- 
lowed by a dichotomy between Ascalabota 
and autarchoglossans (composed of Scinco- 
morpha and Anguimorpha)—many names 
and clade compositions that would be later 
formalized by the first statistical phylogenies 
that used morphological data (see more 
below). Later authors used such essays and 
the morphological characters illustrated in 
them as the basis for cladistic analyses of the 
relationships within groups such as “prim- 
itive” snakes (Rieppel, 1979) and geckos 
(Kluge, 1987). 

A revolutionary classification of snakes 
introduced by Underwood (1967) was based 
on extensive anatomical documentation of a 
nearly exhaustive accounting of both soft 
tissue and osteological characters. In doing 
so, he provided the first systematic overview 
of snake diversity since the 19th century 
works of Cope and Boulenger. He docu- 
mented characteristics of systematic rele- 
vance in extraordinary detail and breadth, 
providing diagnoses or character-based de- 
scriptions for most taxa down to the 
subfamily level, and in many cases to the 
tribe or genus. F urthermore, he explicitly 


recognized the link between such taxonomic 
arrangements, the characters on which they 
are based, and the evolutionary history and 
series of evolutionary transitions they im- 
plied. The enormity of this work and its 
towering importance for nearly all subse- 
quent studies of snake classification are 
difficult to overstate, but we will here only 
note that it formalized the tripartite in- 
fraordinal arrangement of snakes into the 
“lower” snakes, the blindsnakes and allies 
(Scolecophidia) and the boas, pythons, and 
allies (Henophidia), and the “higher” or 
advanced snakes (Caenophidia), such as 
colubrids, dipsadids, elapids, and viperids. 
Despite the extensive amount of research 
in herpetological systematics during the 
middle of the 20th century and the rapid 
expansion of “cladistic” systematics after 
the English publication of Hennig (1966), 
few studies examined any squamate groups 
by explicit quantitative analyses of morpho- 
logical character matrices in the decades 
after Camp (1923). As noted above, only a 
few researchers investigated small groups 
such as “primitive” snakes (Rieppel, 1979) 
and geckos (Kluge, 1987). In contrast, no 
workers attempted an explicit analysis of 
Squamata as a whole for 65 years, until the 
first comprehensive phylogeny of the group 
was estimated by Estes et al. (1988), 
marking the modern period of squamate 
systematics a mere 33 years ago. 


Early Morphological Phylogenies 


The first attempt at a comprehensive 
phylogenetic analysis of broadscale squa- 
mate relationships was conducted by Estes 
et al. (1988). The analysis included 148 
morphological characters and 19 extant 
squamate families as tips (= operational 
taxonomic units) analyzed by maximum 
parsimony—an approach followed by most 
subsequent morphological analyses. Squa- 
mates were diagnosed by several characters, 
including fused premaxillae, fused parietals, 
opisthotics fused to exoccipitals, loss of 
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palatine teeth, procoelous vertebrae, loss of 
posterior trunk (thoracolumbar) intercentra, 
and presence of an anterior coracoid 
emargination. (See Fig. 2 for a brief sample 
of morphological characters that have been 
historically used to diagnose squamates and 
several squamate clades.) Within squa- 
mates, the results of Estes et al. (1988) 
recovered a broad-level structure very 
similar to the highly influential precladistic 
classification of Camp (1923) (Fig. 3A). 
Specifically, the earliest branching event in 
the squamate tree would be that between 
Iguania (diagnosed in part by fleshy, lobular 
tongues; frontals strongly constricted be- 
tween orbits; postfrontals reduced or lost; 
and fingerlike angular process on retroar- 
ticular process) and Scleroglossa (i.e., all 
other squamates, diagnosed in part by scaly, 
forked tongues; prey prehension exclusively 
mediated by jaws rather than tongue; 
postfrontals forked medially; enlarged ce- 
phalic scales; loss of the dorsal process of 
squamosal; alar process of prootic elongat- 
ed; and clavicles strongly angulated, curving 
anteriorly; e.g., Fig. 2, chars. 1, 4, 5, and 11). 

Within scleroglossans, Estes et al. (1988) 
estimated a major dichotomy between 
gekkotans (diagnosed in part by the loss of 
the lacrimals and postorbitals; reduction of 
the jugals and postfrontals; and persistent 
notochordal canal in adults, with amphi- 
coelic vertebrae in some families) and 
autarchoglossans (diagnosed by only three 
synapomorphies: loss of jugal-squamosal 
contact on supratemporal arch; presence 
of dermal bone rugosities; and presence of 
an m. rectus abdominis lateralis muscle). 
Autarchoglossans were subdivided into two 
main lineages: “Scincomorpha and Angui- 
morpha. Scincomorphs included xantusiids, 
lacertids, teiioids (teiids and gymnophthal- 
mids), cordyloids (cordylids and gerrhosaur- 
ids), and scincids, whereas anguimorphs 
comprised anguids, varanids, heloderma- 
tids, and xenosaurids (Fig. 3A). 

In the study by Estes et al. (1988), all 
sampled elongate and limb-reduced or 


legless taxa except pygopodids—snakes, 
amphisbaenians, and dibamids—were 
“clumped” together within anguimorphs. 
However, this was considered a problematic 
result because the group was primarily 
supported by characters associated with 
limb reduction, which were thought to be 
convergent. Therefore, the final and pre- 
ferred phylogenetic tree reported by the 
authors did not indicate the placement of 
those three limb-reduced clades. Pygopo- 
dids, on the other hand, were recovered as 
the sister clade to gekkotans instead of 
clustering with other limb-reduced clades, 
most likely driven by soft tissue characters 
and gecko-specific cranial characters that 
were able to overcome the signal from 
limbless characters. 

A series of subsequent modifications and 
additions to this initial morphological data 
set followed in the years afterwards, includ- 
ing the first studies to add fossil taxa to the 
original all-extant data set (e.g., Wu et al., 
1996; Evans and Barbadillo, 1998; Evans 
and Chure, 1998; Caldwell, 1999). All of 
those studies except Caldwell (1999) found 
a very similar basic configuration related to 
the works of Camp (1923) and Estes et al. 
(1988). Importantly, the latter studies pro- 
vided the first phylogeny-based systematic 
classification of important extinct lineages, 
such as the placement of Borioteiioidea (= 
Polyglyphanodontidae) as a sister clade to 
either Teiioidea or Lacertoidea (Teiioi- 
dea+Lacertidae), and paramacellodids as 
either a sister taxon to or a member of 
Scincoidea (Scincidae+Cordyloidea) (Evans 
and Barbadillo, 1998; Evans and Chure, 
1998). However, the relationships of other 
fossil species varied considerably across 
some of those analyses, such as the place- 
ment of the oldest articulated fossil lizards 
known at the time from the Late Jurassic of 
Germany, including Eichstaettisaurus, Ba- 
varisaurus, and Ardeosaurus (Evans and 
Barbadillo, 1998; Evans and Chure, 1998). 

The largest reassessment of the original 
matrix of Estes et al. (1988) was provided by 
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Figure 3. Competing hypotheses of squamate phylogeny from morphological data. (A) First computer-based morphological 


phylogeny from Estes et al. (1988). (B) Tree from Lee (1998), 


similar to its subsequent expansions (Lee, 2000; 2005a,b; Lee and 


Caldwell, 2000). In this hypothesis, removal or down-weighting characters related to fossoriality/limb reduction is able to separate 
limb-reduced lineages. (C) Summary version of the tree from Conrad (2008), the first to provide intensive sampling of both extant 
and fossil taxa. (D) Summary version of the tree from Gauthier et al. (2012), which included fewer taxa than Conrad (2008) but 


considerably more morphological characters. 


Caldwell (1999), who revised the characters 
from the original study on the basis of their 
interdependence and also removed soft 
tissue characters (which could not be scored 
for fossils), resulting in the exclusion or 
merger of nearly one-third of the original 
characters. New fossil taxa were also includ- 


ed, such as the oldest known fossil snake at 
the time (Dinilysia patagonica) and Mosa- 
sauria (dolichosaurids, aigialosaurids, and 
mosasaurids). The results still estimated 
iguanians as the earliest diverging squamate 
crown clade along with a monophyletic 
Anguimorpha and Scincomorpha. However, 
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gekkotans were found to be the sister 
lineage of scincomorphs instead of autarch- 
oglossans. Importantly, mosasaurs were 
found to be the sister lineage of snakes, 
with Mosasauria+Serpentes representing an 
early-branching lineage in squamate evolu- 
tion. 

After Estes et al. (1988), the first new 
morphological data set aimed at a compre- 
hensive squamate phylogeny was provided 
by Lee (1998). After a considerable reeval- 
uation of previous data sets and inclusion of 
a large number of new observations, this 
matrix totaled 230 characters scored for 22 
taxa. This data set and its subsequent 
modifications and expansions (Lee and 
Caldwell, 2000; Lee, 2005a,b) provided (at 
that time) the most substantial shift from 
the classical framework previously estab- 
lished by Camp (1923) and Estes et al. 
(1988) (Fig. 3B). 

This phylogeny was the first to estimate 
scincomorphs to be paraphyletic: Scincids 
and cordyloids were more closely related to 
anguimorphs (clade Diploglossa) instead of 
forming a clade with lacertoids (lacer- 
tidsHteiioids). This novel clade was diag- 
nosed in part by the retroarticular process 
twisted distally, broad distally and with a 
smooth dorsal surface; osteoderms present 
both dorsally and ventrally on the body and 
over the skull table; and caudal transverse 
processes converging distally. Furthermore, 
Lee (1998) recovered Gekkota and Xantu- 
siidae as sister lineages and within the same 
clade as Dibamidae and Amphisbaenia— 
clade Nyctisaura, a clade name resurrected 
from Cope (1900)—diagnosed in part by the 
absence of a postorbital, an occipital con- 
dyle weakly or strongly bipartite; medial 
closure of the Meckelian canal by the 
dentary; and absence of the angular as a 
discrete element. Finally, this study also 
found Mosasauria to be the sister taxon to 
early-diverging snakes within the clade 
Pythonomorpha (diagnosed by several syn- 
apomorphies, including the quadrate being 
suspended mostly by the supratemporal; 


strongly developed parietal downgrowths; 
stapedial footplate surrounded by the pro- 
otic and opisthotic [but see for the homol- 
ogy of this character Rieppel and Kearney, 
2002]; supraoccipital sutured to parietal; 
mobile contact between angular and sple- 
nial [contributing to an intramandibular 
joint]; recumbent replacement teeth; reduc- 
tion of median premaxillary teeth; and 
presence of zygosphenes and zygantra on 
the vertebral neural arches; e.g., Fig. 2, 
char. 9). 

The resurrection of Pythonomorpha, also 
a name first proposed by Cope ( 1869), for 
grouping mosasaurs and snakes was also 
recovered by other authors around the same 
time (Scanlon, 1996; Caldwell, 1999), in- 
cluding slightly earlier versions of the Lee 
(1998) data set (Lee, 1997a,b). This is 
perhaps the most controversial result among 
morphological phylogenies because it led to 
a series of subsequent debates on the 
marine, terrestrial, or fossorial origin of 
snakes and their placement within the 
squamate Tree of Life (e.g., Rieppel and 
Zaher, 2000; Tchernovy et al., 2000; Conrad, 
2008; Gauthier et al., 2012: see below under 
“The Origin of Snakes”). Those controver- 
sial (although not necessarily unexpected) 
results were subsequently maintained by 
analyzing soft tissue data alone (Lee, 2000), 
and some recent morphological studies (see 
below). It is perhaps surprising that the 
placement of snakes and mosasaurs as sister 
taxa within squamates (both groups never 
having been tested together in a phyloge- 
netic analysis before the late 1990s in the 
aforementioned studies) sparked more con- 
troversy than the proposed paraphyly of 
scincomorphs (a clade that had been 
established for nearly 75 years) and the 
subsequent placement of cordylids and 
scincids with anguimorphs, in direct conflict 
with the traditional classifications of Cam 
(1923) and Estes et al. (1988). This situation 
underscores the effect that forces such as an 
organism's charisma or even its penetration 
into popular culture and religion can have 
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on the broader effect of scientific discovery 
and research agendas, regardless of how 
much they actually contradict previously 
established hypotheses. 


Large-scale Morphological Data Sets 

Two decades after the first quantitative 
phylogenetic analysis of squamate relation- 
ships (Estes et al., 1988), new morphological 
data sets encompassing an increasingly 
larger number of characters and taxa began 
to emerge. These data reflected a growing 
body of evidence during the late 1990s and 
early 2000s regarding the advantages of 
increased taxon and character sampling in 
phylogenetic inference (Graybeal, 1998; 
Wiens, 2004), specifically regarding the 
benefits of including fossil taxa (Wiens, 
2003, 2006). These matrices also represent- 
ed further attempts to “test” the morpho- 
logical signal in squamates (i.e., the strength 
of the morphological data for driving broad- 
level squamate phylogenetic relationships), 
given the contrasting results obtained by the 
first emerging molecular phylogenies, which 
contested all previous morphology-based 
hypotheses of squamate evolution (Town- 
send et al., 2004; Vidal and Hedges, 2005; 
see below). 

The first large-scale morphological data 
set for squamates was that of Conrad (2008), 
including 223 taxa and 363 primarily 
osteological characters. This data set and 
its subsequent expansions and modifica- 
tions, including combined evidence analyses 
(e.g., Wiens et al., 2010: Simdes et al., 
2015a; Conrad, 2017; Pyron, 2017) still 
represent the most extensive taxonomic 
sampling with morphological data available 
today. The overall structure of the tree in all 
studies analyzing the Conrad (2008) matrix 
alone is similar to that of Estes et al. (1988), 
including the early divergence of iguanians 
and Scleroglossa, with a monophyletic 
Autarchoglossa, Scincomorpha, and Angui- 
morpha (Fig. 3C). However, the much 
broader taxonomic sampling of extant and 


fossil taxa had the dual effects of breaking 
long branches and introducing more termi- 
nals with observed character state reversals. 
Thus, considerably fewer characters were 
found to diagnose each major clade unam- 
biguously. For instance, total group Scle- 
roglossa (or Scincogekkonomorpha under 
Conrad’s_ terminology) was diagnosed by 
only five synapomorphies (compared with 
27 by Estes et al., 1988). Similar to previous 
morphology-based phylogenies focused on 
iguanian relationships (Frost and Etheridge, 
1989), Conrad (2008) obtained very low 
resolution and stability among the different 
clades of iguanians, a limiting aspect also 
observed later in Gauthier et al. (2012). 
Snakes, amphisbaenians, dibamids, limb- 
reduced skinks (Feylinidae and Acontidae) 
were clustered in a limb-reduced clade: 
Scincophidia (except pygopodids). Within 
snakes, scolecophidians were the earliest 
evolving clade of snakes, with subsequent 
divergence of extinct snake taxa (e.g., 
Dinilysia, Pachyrhachis, and Wonambi), 
and finally caenophidians+henophidians. 
The overall structure of Conrad’s (2008) 
parsimony tree was later sustained under 
both maximum parsimony and Bayesian 
inference (Wiens et al., 2010). 

The much greater sampling of fossil taxa 
by Conrad (2008) relative to all other 
studies also provided the first assessment 
of some fossil families in a comprehensive 
estimate. For instance, necrosaurids were 
recovered as early-diverging varanoids, 
whereas several fossil species from the Late 
Cretaceous of Mongolia (e.g., Gobiderma, 
Estesia, and Paraderma) were recovered as 
closely related to the Gila monster (helo- 
dermatids) in a new clade named Monster- 
sauria. Importantly, Mosasauria was 
recovered as deeply nested within varanoids 
and only distantly related to snakes, in sharp 
contrast to most of the work done during 
the previous decade that established a sister 
group relationship between snakes and 
mosasaurs (Caldwell and Lee, 1997; Lee, 
1998: Caldwell, 1999; Lee and Caldwell, 
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2000; Lee, 2005b; but see Rieppel and 
Zaher, 2000). 

The next step towards the expansion of 
character and taxon sampling for squamate 
phylogeny was provided by Gauthier et al. 
(2012), with 192 taxa and 610 morphological 
characters, including highly detailed char- 
acter descriptions and illustrations. Similar 
to Conrad (2008), despite a considerable 
increase in the number of taxa and charac- 
ters, the overall tree structure of extant 
squamate families under maximum parsi- 
mony reinforced the framework first intro- 
duced by Camp (1923) and Estes et al. 
(1988) (but different from Lee [1998] and 
its subsequent expansions) and contested 
results obtained from molecular data sets 
(Fig. 3D). However, Bayesian inference 
revealed a significant lack of resolution by 
finding a large polytomy among the major 
iguanian and scleroglossan clades. As with 
nearly all previous morphology-based squa- 
mate phylogenies, most limb-reduced taxa 
(except pygopodids) formed a clade that was 
inclusive of snakes, amphisbaenians, and 
dibamids, plus Anniella (a limb-reduced 
anguid) and the fossil taxon Sineoamphis- 
baena. This assemblage was called a “fos- 
sorial group,” although the vast majority of 
snake species are not fossorial (see Wallach 
et al., 2014). 

Gauthier et al. (2012) included a much 
larger sampling of snake taxa than any 
previous studies and found contrasting 
results concerning early snake phylogeny. 
Whereas Cretaceous fossil snakes such as 
Dinilysia and Najash were estimated as the 
earliest evolving snakes under maximum 
parsimony, scolecophidians were found as 
the earliest evolving snakes under Bayesian 
inference. Also in contrast to previous 
studies, Gauthier et al. (2012) found new 
and controversial placements for some fossil 
lizard lineages, most notably borioteiioids 
and mosasaurians. The two groups were 
estimated as non-scleroglossan squamates, 
suggesting that they are the earliest diverg- 
ing squamate lineages after the initial 


divergence of Iguania, and thus successive 
sister lineages to Scleroglossa. The subse- 
quent expansion of their data set by Reeder 
et al. (2015) with the addition of 81 soft 
tissue characters (totaling 691 morphologi- 
cal characters), yielded similar results when 
analyzed alone under all tested optimality 
criteria (maximum parsimony, maximum 
likelihood, and Bayesian inference). 


Enduring Controversies 


After nearly 100 years of attempts to 
establish relationships among the major 
squamate lineages, dozens of different 
analyses and data sets have yielded a few 
consistent but occasionally opposing para- 
digms with the use of morphological data. 
First, most studies agree that many squa- 
mate groups long recognized by prephylo- 
genetic workers (going back to Linnaeus or 
earlier) are monophyletic, including snakes, 
geckos, skinks, and several fossil lineages, 
including Mosasauridae, Borioteiioidea, 
Paramacellodidae, Madtsoiidae, and Glyp- 
tosaurinae. A primary challenge, however, 
has been identifying how those different 
groups are related to each other. Another 
common point of agreement is in the 
placement of Iguania as the earliest diverg- 
ing group of crown squamates. There is also 
widespread agreement on a close relation- 
ship between anguids, helodermatids, xen- 
osaurids, and varanoids, comprising 
Anguimorpha. Most studies also concur on 
a close relationship between paramacello- 
dids and scincoids (either within scincoids 
or the sister clade), as well as the placement 
of Mosasauria within anguimorphans, re- 
gardless of whether snakes were found as 
the sister group to Mosasauria or elsewhere. 

The most considerable sources of dis- 
agreement between morphological hypoth- 
eses stem from two major groups of 
competing data sets and analyses: one with 
the basic tree structure first introduced by 
Camp (1923) and maintained by Estes et al. 
(1988), Conrad (2008), and Gauthier et al. 
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(2012) and the one introduced by Lee 
(1998) and further promoted by Lee and 
Caldwell (2000), among others. Whereas 
the former is in agreement on the mono- 
phyly of Scincomorpha and the placement 
of snakes in the “fossorial” clade, including 
dibamids and amphisbaenians, the latter 
suggests the placement of scincoids in a 
clade with Anguimorpha and the placement 
of snakes as the sister lineage to Mosasauria. 
Another and less-recognized source of 
discordance between morphological data 
sets is the conflict regarding the placement 
of geckos. Whereas most morphological 
phylogenies recovered gekkotans as the 
sister clade to autarchoglossans, diverging 
early in squamate evolutionary history 
subsequent to iguanians (e.g., Estes et al., 
1988; Conrad, 2008; Gauthier et al., 2012), 
other morphological data sets recovered 
gekkotans in a clade with Mosasauria and 
snakes (Caldwell, 1999; Simoes et al., 
2017b), close to borioteiioids (Lee, 2005a), 
or possibly forming a clade with borio- 
teiioids or amphisbaenians, dibamids, and 
xantusiids, or a combination of these groups 
(Lee, 2000, 2005a,b). 

Similarly, the controversy over the place- 
ment of snakes as either the sister lineage of 
Mosasauria or in a clade with other limb- 
reduced taxa sparked a wave of debate over 
the origin of snakes, including marine (Lee, 
1997a, 1998, 2005a,b), terrestrial (Vidal and 
Hedges, 2004), or fossorial (e.g., Rieppel 
and Zaher, 2000; Gauthier et al., 2012: 
Hsiang et al., 2015) hypotheses for the 
ecological driver of early snake evolution. 
The “snake origins” debate further boosted 
another discussion on the reliability of 
morphological characters in the face of 
widespread apparent homoplasy, particular- 
ly those associated with limb reduction in 
squamates (Lee, 1998, 2005a,b). Characters 
highly correlated with limb reduction, 
miniaturization of the skull, and fossorial 
adaptations seem invariably to “attract” 
limb-reduced taxa to cluster in morpholog- 
ical data sets, unless each limb-reduced 


taxon is included and analyzed separately 
from other limb-reduced clades or character 
weighting is applied (Lee, 1998, 2005b). 
Yet, few if any researchers actually “be- 
lieve” that these taxa form a natural group, 
particularly when it includes species such as 
legless skinks, which are almost indisputably 
highly derived and deeply nested in Scinci- 
dae, for instance (see Wiens and Lambert 
[2014] for more discussion of this point). 
Subsequent studies eventually called into 
question the entire conceptual and meth- 
odological approach taken during the con- 
struction of nearly all previous squamate 
morphological data sets. Upon a reassess- 
ment of more than 1,000 morphological 
characters used by earlier authors, Simdes 
et al. (2017b) found that ~35%—-45% of all 
characters used in previous analyses violat- 
ed basic principles of morphological char- 
acter construction, including not accounting 
for logical and biological dependencies 
among characters (Wilkinson, 1995; Strong 
and Lipscomb, 1999), not keeping multiple 
exclusivity among states of a single character 
(Sereno, 2007; Brazeau, 2011), violating 
principles of character coding that avoid 
biases introduced by inapplicable or missing 
data (Strong and Lipscomb, 1999; Brazeau, 
2011), and not treating continuous variables 
as such by introducing arbitrarily delimited 
character discretization (Rae, 1998: Golob- 
off et al., 2006). Upon reevaluating the 
largest available squamate morphological 
data sets at the time (Conrad, 2008; 
Gauthier et al., 2012) and recoding or 
removing problematic characters, Simdes 
et al. (2017b) found results in substantial 
contrast to the original studies. A reanalysis 
of the matrix of Conrad (2008) resulted in a 
tree with strong similarities to that of Lee 
(1998) and subsequent expansions thereof, 
by estimating scincids and cordyloids as 
more closely related to anguimorphs (Dip- 
loglossa). In a similar reanalysis of Gauthier 
et al. (2012), Mosasauria moved from an 
early-diverging position in the squamate 
tree (an outgroup to Scleroglossa) to form 
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a clade along with the limb-reduced squa- 
mates, including snakes. 

Following stricter criteria for morpholog- 
ical character construction and an expanded 
taxonomic sampling of rhynchocephalians, 
early lepidosaurs, and other diapsid lineag- 
es, Simoes et al. (2018) built a new 
morphological data matrix inclusive of the 
major stem and crown extinct and extant 
lineages of squamates. The result was the 
first ever agreement of the morphological 
signal with the molecular signal (see below) 
for several key nodes in the squamate Tree 
of Life (i.e., gekkotans at the base of the tree 
and iguanians within a monophyletic Tox- 
icofera; Fig. 4). Numerous fossil lineages 
were included, with borioteiioids recovered 
in their classical placement as closely 
related to teiids (specifically, the sister 
group to Teiioidea) and paramacellodids as 
the sister group to skinks. Mosasaurians 
were found to be the sister taxon to snakes 
within Toxicofera, and several previously 
overlooked fossil lepidosaurs were found to 
represent some of the oldest known squa- 
mate fossils, expanding the fossil record of 
the group to the Middle Triassic. This result 
offers great promise for unifying the phylo- 
genetic signal found in different ontological 
partitions of squamate biology (see Losos et 
al., 2012), perhaps the final remaining 
frontier for understanding the squamate 
Tree of Life. 


MOLECULAR HYPOTHESES 
History of the Molecular Hypothesis 


Comprehensive molecular analyses of 
squamate phylogeny were relatively late in 
arrival, only beginning in earnest several 
years into the 21st century. In the 1990s, 
however, a few studies began to suggest that 
molecular data may exhibit substantially 
different phylogenetic signals than tradi- 
tional morphological hypotheses (see 
above). In a study of mitochondrial transfer 
RNAs (tRNAs), Kumazawa and Nishida 
(1995) sampled a gecko, a skink, and an 
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Figure 4. Morphology-only phylogenetic tree from Simoes et 
al. (2018). A new composition of early squamates is recovered, 
with gekkotans as the first branching crown clade of 
squamates and with iguanians forming a clade with 
anguimorphans, snakes, and mosasaurians. 


iguanian, and found that either the skink or 
the gecko represented the earliest diver- 
gence in squamates in various analyses, but 
never the iguanian. In a data set sampling 
the mitochondrial gene ND4 and associated 
tRNAs for several snakes, an anguimorph, a 
gecko, two iguanians, and a skink, Forstner 
et al. (1995) found that the skink was the 
sister lineage to the remaining taxa. A 
reanalysis of this data set by Macey and 
Verma (1997: 277) confirmed support for 
this placement, noting that “no topology is 
consistent with the squamate phylogenetic 
hypotheses derived from morphological 
data.” Finally, Gorr et al. (1998: 481) noted 
that hemoglobin sequences did not place 
Serpentes closer to Anguimorpha than 
Iguania as hypothesized by Estes et al., 
(1988) and others, stating “neither single 
chain nor any of the of tandem 
alignment trees ... are consistent with this 
hypothesis.” 
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Starting with Saint et al. (1998), a series 
of studies (e.g., Harris et al., 1999: Harris, 
2003) began to assemble increasingly well- 


sampled data sets comprising fragments of 


the single-copy, intron-free nuclear protein— 


coding gene CMOS for a broad array of 


Squamate species. In essentially all analyses, 
Iguania, Anguimorpha, and Serpentes 
formed a clade, and Teiioidea was the 
earliest diverging squamate lineage. As 
Harris (2003: 540) noted, “such a relation- 
ship is completely anomalous to estimations 
of relationships based on morphology or the 
fossil record.” Vidal and Hedges (2004) 
subsequently combined parts of these data 
sets with sequence fragments from a second 
protein-coding nuclear gene (RAG1), sam- 
pling nearly all squamate families recog- 
nized at the time. Node support varied, but 
their analyses estimated a topology with 
Dibamidae as the earliest diverging squa- 
mate lineage, followed by geckos; skinks and 
relatives; teiids, lacertids, and amphisbae- 
nians; and a clade comprising Anguimorpha, 
Iguania, and Serpentes (Fig. 5). Despite the 
extraordinary incongruence with previous 
morphological hypotheses, Vidal and Hedg- 
es (2004) were primarily concerned with 
snake origins and thus did not offer much 
additional commentary on the implications 
of their results for understanding higher 
level squamate relationships. 

The landmark study of Townsend et al. 
(2004) was the first multilocus data set 
(CMOS+RAGI1) to address explicitly higher 
level squamate relationships with relatively 
comprehensive taxon sampling. They cor- 
roborated Vidal and Hedges (2004) in 
cementing the “molecular” hypothesis, in 
which dibamids or geckos are the earliest 
diverging lineages, followed by skinks and 
relatives; teiids, lacertids, and amphisbae- 
nians; and, again, a clade comprising 
Anguimorpha, Iguania, and Serpentes (Fig. 
5). They noted “this unconventional rooting 
does not seem to be due to long-branch 
attraction, base composition biases among 
taxa, or convergence caused by similar 
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Figure 5. The phylogenetic tree of squamates according to 
most molecular data sets. Gekkotans and dibamids form the 
first radiation of crown clade squamates, with the placement of 
dibamids still under relative debate. 


selective forces acting on nonsister taxa” 
(Townsend et al., 2004: 735) and suggested 
that the morphological topology was likely 
due to homoplasy. These results were 
quickly corroborated by Vidal and Hedges 
(2005), sampling nine nuclear protein— 
coding genes (albeit for fewer taxa), who 
introduced a tree-based nomenclature for 
these groups. In particular, Scinciformata 
comprises Scincidae, Cordylidae, and Xan- 
tusiidae, whereas Laterata comprises 
teioids, lacertids, and amphisbaenians. Fi- 
nally, the name Toxicofera was applied to 
the clade of Anguimorpha, Iguania, and 
Serpentes. Like Townsend et al., (2004) 
Vidal and Hedges (2005) also suggested that 
phenotypic homoplasy was the most likely 
explanation for the unorthodox results 
compared to the morphological hypothesis. 

A flood of subsequent research yielded 
highly similar results across a variety of data 
sets, including full mitochondrial genomes 
(Kumazawa, 2007), dozens of protein-cod- 
ing nuclear genes (Wiens et al., 2012), 
thousands of species sampled for a smaller 
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set of mitochondrial and nuclear genes 
(Pyron et al., 2013), thousands of ultra- 
conserved element (UCE) loci (Streicher 
and Wiens, 2017), thousands of orthologous 
loci generated from RNA transcriptomes 
(Irisarri et al., 2017), hundreds of anchored 
hybrid-enrichment (AHE) loci (Burbrink et 
al., 2020), and a combination of UCEs, 
AHE, and nuclear protein-coding loci 
(Singhal et al., 2021). This drastic overhaul 
of traditional, primarily morphology-based 
classification mirrors similar situations that 
occurred in other groups, such as lissam- 
phibians (see Frost et al., 2006). Thus, the 
“molecular” hypothesis of squamate rela- 
tionships is now broadly accepted by most 
systematists. Questions remain, however, 
regarding the incongruence with existing 
morphological data sets and the many 
poorly supported branches observed across 
the myriad molecular studies referenced 
above. 


Molecular Data Sets and Uncertainty 


As noted, consistency across molecular 
data sets is remarkable in terms of support 
for the “molecular” hypothesis of a nested 
Iguania (Fig. 5). The signal for the molec- 
ular hypothesis is not marker specific but is 
found across both the nuclear and mito- 
chondrial genomes in both coding and 
noncoding regions. However, support is 
not universal across these loci. For example, 
seemingly artifactual results (paraphyletic 
Iguania and early-branching Serpentes) 
from some whole-mitochondria alignments 
are apparently the result of long-branch 
attraction and an ancient episode of con- 
vergent adaptive evolution (Castoe et al., 
2009). Similarly, the BovB LINE family of 
transposable elements has experienced ex- 
traordinary levels of ectoparasite-mediated 
horizontal transfer between mammals and 
squamates (see Pasquesi et al., 2018, and 
references therein), and its phylogeny is 
essentially uninterpretable for understand- 
ing squamate relationships. 


Almost as remarkable as the consistency 
across studies for the molecular hypothesis 
is their consistency in reflecting the same 
regions of uncertainty in the tree, with a set 
of poorly supported branches being com- 
mon to most analyses (see Burbrink et al., 
2020; Singhal et al., 2021, and references 
therein). Along the backbone of the tree, 
two major nodes merit discussion. The first 
is the most recent common ancestor 
(MRCA) of Squamata, where the earliest 
branch is either Gekkota, Dibamia, or both, 
forming a clade that is the sister lineage of 
remaining squamates. The second is the 
resolution of Toxicofera: whether Serpentes 
is the sister lineage to Iguania, Anguimor- 
pha, or both. 

Within higher level clades, four major 
groups stand out. The first are pleurodont 
iguanians, where relationships among the 
numerous families are famously intractable 
(see Townsend et al., 2011). The second is 
monophyly of Amphisbaenia: whether Lac- 
ertidae is more closely related than Rhine- 
uridae to the remainder of amphisbaenians 
(Wiens et al., 2010; Streicher and Wiens, 
2017). The third is monophyly of Scoleco- 
phidia (blindsnakes), where Anomalepidi- 
dae is sometimes the sister lineage of all 
other snakes exclusive of Typhlopi- 
dae+Leptotyphlopidae, or of all other 
snakes (e.g., Burbrink et al., 2020). Fourth, 
within Caenophidia, relationships among 
the numerous families of booid and colub- 
roid snakes often receive strong support for 
vastly different topologies (see Pyron et al., 
2014; Reynolds et al., 2014; Streicher et al., 
2016, and references therein) mirroring the 
situation in pleurodont iguanians. 

Most of these instances (except Amphis- 
baenia) appear to stem from the well-known 
problem of rapid evolutionary radiations 
and long descendant branches. These sce- 
narios can be difficult if not impossible to 
resolve because of a lack of informative sites 
accumulating during brief periods of diver- 
gence, as well as frequent homoplasy arising 
along deep stem lineages (Rokas and Car- 
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roll, 2006). Phylogenomic data sets alone 
may not be enough to overcome these 
various sources of nonphylogenetic signal 
and other artifacts, which may stem from 
nonorthology, model error for gene trees or 
species trees, taxon sampling, and so on 
(Philippe et al., 2011). In contrast, non- 
monophyly of Amphisbaenia is apparently 
due to long-branch attraction and is seem- 
ingly “corrected” by increased taxon sam- 
pling in molecular data sets (see below and 
Pyron et al., 2013). 


Advancements and Controversies 


Recent attempts to understand the dis- 
tribution and variation of phylogenetic 
signal in phylogenomic data sets for squa- 
mates have corroborated these hypothe- 
sized explanations for topological 
discordance and low support but have 
unfortunately offered limited resolution, 
although with some promising potential 
avenues for future investigation. McMahan 
et al. (2015) suggested that since the branch 
leading to the sole outgroup Sphenodon was 
so long, the signal for the root of Squamata 
was most likely obliterated. Thus, the 
rooting of the ingroup may essentially be 
random. They then claimed to show that 
rerooting the molecular tree with Iguania is 
supported by more molecular synapomor- 
phies than when the tree is rooted with 
Teiioidea or Dibamidae. 

Harrington et al. (2016) showed that this 
procedure erroneously conflated branch 
support with branch length. The number 
of molecular synapomorphies counted along 
a branch is an estimate of the length of that 
branch in substitutions, not a measure of 
character support. Thus, McMahan et al. 
(2015) were in effect arguing that the 
longest ingroup branch of any clade should 
be the root. Furthermore, their counts were 
obscured by collapsing ingroup clades to a 
single basal dichotomy and by considering 
only two possible roots. When the full tree is 
analyzed under their procedure, both 


Teiioidea and Serpentes represent “longer” 
root branches than Iguania. Additionally, no 
known gene tree places Iguania at the root, 
and none of the dozens of published 
molecular phylogenies have obtained this 
result. If the Sphenodon rooting was “ran- 
dom,” one would presumably expect occa- 
sional molecular support for a basal Iguania 
or a non-Gekkota/Dibamia root in at least 
some trees. 

Another attempt to explain the differenc- 
es between (most of) the morphological and 
molecular phylogenetic data sets (Koch and 
Gauthier, 2018) raised similar objections, 
claiming to show that the nested placement 
of Iguania was due to systematic bias and 
noise in molecular data—in essence, long- 
branch attraction of Iguania and Serpentes. 
They purported to demonstrate this by 
removing “fast-evolving” sites, which col- 
lapsed support for Toxicofera and eventually 
all backbone nodes as more characters were 
removed. Subsequently, the authors showed 
that individual genes did not decisively 
support Toxicofera with significance tests 
of topology. They suggested that because 
85% of 46 genes from Reeder et al. (2015) 
included more than one topological alter- 
native in the confidence set, Toxicofera was 
not supported decisively by the molecular 
data. Finally, they offered an analysis of 
phylogenetic informativeness (PI) of the 46 
genes from Reeder et al. (2015). These 
showed younger peaks in PI than the ages of 
the root branches of interest, concluding 
that these data sets thus had inadequate 
power to resolve early divergences. 

Owing to the recency of this criticism, to 
our knowledge it has so far not been 
critically evaluated in the literature. How- 
ever, some immediate limitations are ap- 
parent in their conclusions, including: 


1) The authors did not provide any direct 
evidence that “fast-evolving” sites were 
artifactually biased by any specific mech- 
anism, beyond showing that snakes and 
iguanians have higher rates of molecular 
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evolution and a slight Adenine-Thymine 
(AT) bias. Thus, they have shown only 
that removing the most potentially in- 
formative sites reduces phylogenetic 
resolution, which is a truism that does 
not bear on the placement of Iguania. 
Importantly, it has been recognized that 
the branch leading to snakes also exhibits 
exceptionally high rates of morphological 
evolution (Watanabe et al., 2019; Sim6es 
et al., 2020b) which are in fact propor- 
tionately faster than molecular rates of 
evolution (Sim6es et al., 2020b). 

The lack of decisiveness for individual 
genes is a theoretical expectation for 
short internodes because of incomplete 
lineage sorting that is typically resolved 
by the use of species tree methods 
(Edwards, 2009) and thus does not 
reflect on the inclusive phylogenetic 
signal of the constituent genes. In 
contrast, well-sampled data sets of sto- 
chastically discordant gene trees are 
expected to converge on the correct 
species tree with high confidence, even 
when a low proportion of individual 
genes matches the true topology (Ed- 
wards et al., 2007; Edwards, 2009). 
Although the root nodes of Squamata do 
not exactly match the peak of PI in terms 
of timing, they are still close to it, with 
drastically higher PI than later in the 
tree. Under the authors’ logic, it should 
thus be difficult to estimate recent 
divergences in Squamata with these 
data; yet, as noted earlier, morphological 
and molecular data are highly congruent 
for many younger nodes and estimated 
relationships. It is also unclear, under 
the authors’ formulation, how it would 
ever be possible to estimate more 
ancient divergences with such data 
(e.g., Hugall et al., 2007), in that they 
come before the peak of PI. 

Although the authors briefly mention the 
4.178 UCE data set of Streicher and 
Wiens (2017), they only tested the signal 
of the 46 loci of Reeder et al. (2015). 


Thus, they did not address how the 
stochastic factors that they claim drive 
the putatively erroneous estimate of 
Toxicofera could be so consistently and 
widely distributed across both the nucle- 
ar and mitochondrial genome and thou- 
sands of coding and noncoding loci, 
without even occasional support for a 
basal Iguania being found among thou- 
sands of independent gene trees. 


More recently, two studies took a sys- 
tematic approach to examining locus-level 
variation in the phylogenetic signal in 
squamate genome-scale data sets, one using 
394 AHE loci (Burbrink et al., 2020) and a 
combined probe set (Singhal et al., 2021) 
consisting of 38 nuclear protein-coding loci 
from Wiens et al. (2012), 5,052 UCEs from 
Streicher and Wiens (2017), and 372 AHE 
loci from Singhal et al. (2017). Both were 
consistent in recovering a reassuring set of 
explanatory variables for discordance among 
loci. Generally speaking, longer loci with 
more informative sites are more congruent 
in supporting the molecular hypothesis and 
Toxicofera. Many loci that do not support 
these nodes instead exhibit incomplete 
lineage sorting, which is resolved by the 
use of coalescent species tree methods and 
subsequently favors Toxicofera. Finally, 
some discordance is due to short internodes 
in the species tree, for which few informa- 
tive sites have accumulated in any locus. 
Additionally, some gene tree error is driven 
by noise arising from a nonphylogenetic 
signal such as loci with few informative 
sites, model violations such as rate hetero- 
geneity or saturation, and alignment issues 
such as low sequencing or assembly quality 
or incorrect ortholog identification (see 
Phillippe et al. 2011). However, these low- 
quality loci do not support any alternative 
topology strongly, and their removal only 
increases support for the molecular hypoth- 
esis and Toxicofera, while no locus supports 
Scleroglossa or a basal Iguania. 
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Concomitantly, both Burbrink et. al. 
(2020) and Singhal et al. (2021) reported 
continuing uncertainty regarding numerous 
nodes. Chief among these is the placement 
of Dibamia, which represents a very short 
internode subtending a very long branch; 
conditions that make resolution extremely 
difficult (Philippe et al., 2011) and may 
continue to be intractable regardless of the 
number of loci sampled until other types of 
data are examined (see below), Similarly, 
the resolution of Toxicofera remains inde- 
cisive, although most recent analyses are 
congruent in estimating an Iguania- 
Anguimorpha clade. The same is true of 
scolecophidian monophyly, where many 
recent data sets are congruent in recovering 
Anomalepididae as the sister lineage to all 
other snakes exclusive of Typhlopoidea. 
Like Dibamia, both Toxicofera and Scole- 
cophidia represent short internodes with 
long descendent branches that may require 
alternative data sources to resolve conclu- 
sively. 

Another contentious point is monophyly 
of Amphisbaenia. Burbrink et al. (2020) did 
not sample Rhineuridae (likely the earliest 
branch in the group), but both Streicher 
and Wiens (2017) and Singhal et al. (2021) 
included it and estimated weak support for 
amphisbaenian paraphyly. This is an appar- 
ent case of long-branch attraction that has 
been putatively solved through increased 
taxon sampling in other molecular data sets 
(see below). Taxon sampling may be a factor 
in many of the examples discussed here, 
because most genome-scale data sets thus 
far have only sampled dozens or a few 
hundred representatives from the ~11,000 
species total diversity of Squamata. 

Finally, the problematic radiations of 
booid snakes and pleurodont iguanians were 
shown to be in the anomaly zone (Singhal et 
al., 2021), where short internodes predom- 
inantly produce discordant gene trees and 
cause nearly insurmountable biases for 
estimating the species tree (Degnan and 
Rosenberg, 2006). Although sampling a 


greater number of informative loci should 
eventually allow resolution in these scenar- 
ios (Linkem et al., 2016), that level has 
clearly not been reached, even with thou- 
sands of markers in squamates. As few as 
three short internodes in rapid succession 
can produce as many as 45 anomalous gene 
trees (Rosenberg and Tao, 2008), and 
Burbrink et al. (2020) and Singhal et al. 
(2021) find at least four and nine short 
internodes in Pleurodonta, respectively. 
Thus, whole-genome data or other charac- 
ter types (e.g., synteny mapping; see below) 
may be needed to provide any conclusive 
resolution for these branches. Neither 
Burbrink et al. (2020) nor Singhal et al. 
(2021) sampled enough colubroid or elapoid 
snake families to make robust inferences 
regarding those historically difficult groups 
(see Pyron et al., 2014), but similar pro- 
cesses seem likely to be at play, with 
numerous short branches during rapid 
Cenozoic radiations in advanced snakes 
(Pyron et al., 2011; Zaher et al., 2019). 


COMBINED EVIDENCE OF SQUAMATE 
RELATIONSHIPS 


Early Searches for Congruence 


The phylogeny of Estes et al. (1988) 
cemented the “morphological” hypothesis 
in the minds of most herpetological system- 
atists as the best quantitative estimate of 
squamate relationships. Concomitantly, as 
molecular data increasingly gained ascen- 
dancy in the field of systematics during the 
late 1990s and early 2000s, such data sets 
provided spectacular new resolutions of 
phylogeny and taxonomy for numerous 
squamate groups. Studies were published 
on dozens of families and genera, with 
examples such as Gymnophthalmidae (Pel- 
legrino et al., 2001), rat snakes in the genus 
Elaphe sensu lato (Utiger et al., 2002), 
Scincidae (Whiting et al., 2003), Iguania 
(Sites et al., 1996; Frost et al., 2001), and 
dwarf boas (Wilcox et al., 2002), represent- 
ing merely a random sample of papers from 
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this period with >100 citations on Google 
Scholar as of February 2021. Anticipation in 
the field for a widely sampled phylogeny of 
Squamata thus began to increase rapidly, 
holding the potential promise of merging 
morphological and molecular data matrices 
to provide a new, unified view of lepidosaur 
evolution. 

At the same time, molecular trees such as 
those of Saint et al. (1998) began to 
accumulate in the literature, and rumors 
of preliminary results from several labs that 
would result in publications such as Town- 
send et al. (2004) and Vidal and Hedges 
(2004) began to circulate among workers in 
the field. Thus, many herpetologists gradu- 
ally became cognizant that the “molecular” 
tree might end up being significantly 
different from the “morphological” tree. 
With this dawning awareness in the back- 
ground, the United States National Science 
Foundation (NSF) awarded approximately 
US$2.4 million in 2004 to an international 
team of researchers under the “Assembling 
the Tree of Life” program to investigate 
squamate origins with molecular and mor- 
phological data and estimate a total evi- 
dence phylogeny. Titled “The Deep Scaly 
Project: Resolving Squamate Phylogeny 
Using Genomic and Morphological Ap- 
proaches,” this project aimed at a grand 
synthesis. As the official NSF abstract stated 
(DEB-0334966): 


Many critical questions in squamate 
evolution remain unresolved, such as 
identification of the most primitive line- 
age of squamates, the origin of snakes, 
and the relationships of venomous snakes 
to other snake lineages. An international 
team of eight investigators from diverse 
institutions (Brigham Young University, 
Field Museum of Natural History, San 
Diego State University, State University 
of New York-Stony Brook, University of 
Adelaide, University of Texas-Austin, and 
Yale University) will collaborate to re- 
solve squamate relationships. Anatomical 


data from living and fossil forms will be 
combined with DNA sequences from 50 
genes for 142 representative squamate 
species. Anatomical data will be obtained 
using traditional methods and new high 
resolution X-ray scanning techniques. 
DNA data will be generated by incorpo- 
rating new tools and databases from 
recent vertebrate genome projects. Com- 
puter modeling will be used to deter- 
mine how data from molecular and fossil 
studies can best be combined to recon- 
struct evolutionary trees. 


As the numerous publications referenced 
above appeared between 2004 and 2010, it 
became increasingly apparent that an illu- 
minating evolutionary consensus of the 
various data sets would not be forthcoming. 
After 8 years of attempting to find a 
congruent explanation linking the molecular 
and morphological data sets, Daza (2014: 
341) notes that: “... the two teams diverged 
and published independently, and each 
group developed a different phylogeny 
based on their own independent analyses 
of a data set that was composed of the 
largest data partitions for squamate reptiles 
ever produced (Gauthier et al., 2012; Wiens 
etal, 2012)> 

In many ways, this fundamental diver- 
gence is still among the primary questions 
for researchers in squamate systematics. As 
Losos et al. (2012) asked: Who speaks with a 
forked tongue? Which data set is “lying” 
about the true order of branching events? 
They note (Losos et al., 2012: 1429): “When 
two phylogenies are fundamentally discor- 
dant, at least one data set must be 
misleading. ... We are left with a conun- 
drum. The molecular data imply an aston- 
ishing pattern of morphological homoplasy 
and suggest very limited knowledge of the 
functional link between structures and 
lifestyle; if convergence is so pervasive, 
what faith can we have in the placement 
of fossil taxa for which no molecular data are 
available? Conversely, morphology implies a 
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pattern of molecular evolution that has yet 
to be explained.” 

It is from this basis, a fundamentally 
puzzling and apparently intractable clash of 
phylogenetic signals from what should be 
epistemologically equivalent data sets, that 
all subsequent attempts at integration and 
consilience have proceeded. These include 
early attempts published before the main 
Deep Scaly data sets (Lee, 2005a, 2009), 
continuations of the Deep Scaly project 
(e.g., Reeder et al. 2015), comparisons of 
multiple existing morphological data sets 
including Conrad (2008) and Gauthier et al. 
(2012) by Pyron (2017), and de novo 
construction of morphological matrices in 
an attempt to understand and eliminate bias 
in character ontology and coding (SimGes et 
al., 2017b, 2018). With this recent historical 
background in mind, we now examine some 
of these results in more detail. 


Results, Successes, and Limitations 


To our knowledge, the first study to 
combine morphological and molecular data 
as independent partitions to assess broad- 
scale relationships in squamates was that of 
Lee (2005a). In that study, the molecular 
data set from Vidal and Hedges (2005) was 
combined with an updated version of a 
previous morphological data set (Lee and 
Caldwell, 2000). Given the relatively small 
size of the molecular data set (especially in 
terms of informative sites), the combined 
evidence tree analyzed under maximum 
parsimony greatly resembled the framework 
obtained by analyzing the morphological 
data set of Lee (1998) and its subsequent 
expansions alone. Specifically, the earliest 
split in the squamate tree represented by 
that between iguanians and scleroglossans, a 
monophyletic Nyctisaura (Xantusii- 
dae+Gekkota+Dibamidae+A mphisbaenia), 
and Diploglossa (Scincoidea+Angui- 
morpha), in which early mosasaurians rep- 
resent the sister taxa to Serpentes. Some 
important novelties however included the 


placement of borioteiioids and lacertoids as 
the earliest diverging branches within Scle- 
roglossa, before the dichotomy between 
Nyctisaura and Diploglossa. 

Despite the initial dominance of the 
morphological signal over the molecular 
signal, this scenario soon reversed with the 
incorporation of larger molecular data sets 
into combined evidence studies of squamate 
phylogeny. After combining nine nuclear 
and one mitochondrial loci (totaling 8,086 
bp) from both Townsend et al. (2004) and 
Vidal and Hedges (2005), Lee (2009) found 
a result that was largely congruent with 
previous molecular phylogenies of squa- 
mates. Ever-growing molecular data sets 
increasingly tended to override most signal 
available from morphological partitions. The 
molecular component of combined evi- 
dence data sets increased to 22 nuclear loci 
(totaling 15,794 bp) analyzed with the 363 
morphological characters from Conrad 
(2008) by Wiens et al., (2010), and then to 
46 loci (45 nuclear and 1 mitochondrial, 
totaling 35,673 bp) with 691 morphological 
characters (Reeder et al., 2015). These 
studies consistently found the estimated 
trees to reflect the overall structure of the 
molecular signal, albeit with widely variable 
placements of the fossil species and extinct 
lineages. 

Similar effects were observed in subse- 
quent studies implementing smaller mo- 
lecular data sets, including Pyron (2017) 
with six nuclear loci (totaling 8,710 bp) and 
a morphological partition of either 363 or 
610 characters from Conrad (2008) and 
Gauthier et al. (2012), and Simées et al. 
(2018) with 16 loci (13 nuclear and 3 
mitochondrial, totaling 11,532 bp) and 347 
newly assembled morphological characters 
(Fig. 6). However, we note that the number 
of changes from the morphological to the 
combined evidence tree in Simdées et al. 
(2018) is somewhat reduced because the 
morphological tree from the latter already 
recovers relationships more similar to the 
molecular tree (Fig. 4). Furthermore, 
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Figure 6. One recent example of the phylogenetic tree of squamates inferred from combined morphological and molecular data 


from Simoes et al. (2018)— 


branch lengths are proportional to chronological time. Relationships among the major squamate clades 


is strongly influenced by the molecular signal (Fig. 3), and similar topologies are obtained by Wiens et al. (2010), Reeder et al. 
(2015), and Pyron (2017), with the majority of differences lying on the placement of some fossil lineages. Abbreviations: Gek, 


Gekkota; Scin, Scincoidea; Serp, total group Serpentes. 


Wiens et al. (2010) and Pyron (2017) noted 
cases in which the addition of morpholog- 
ical data altered the molecular hypotheses, 
as well, indicating that even relatively small 
morphological partitions can still exert a 
decisive influence on larger molecular 
matrices. 


Generally, the overwhelming signal from 
molecular data in nearly all combined 
evidence analyses has an indirect effect 
even on the placement of extinct lineages. 
The placement of fossil taxa can only be 
directly influenced by morphological data. 
However, molecular data influence the 
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relationships between extant lineages for 
which morphological data are also available. 
Thus, molecular data can thereby indirectly 
influence the placement of fossil taxa by 
altering the likelihood of various arrange- 
ments of branches relative to the morpho- 
logical matrix. The reverse is also true in 
theory, although such effects are typically 
observed to a far lesser extent in most 
empirical analyses (see Wiens et al.. 2010: 
Pyron, 2017). 

For instance, Huehuecuetzpalli_ mixtecus 
(Early Cretaceous, Mexico) is commonly 
estimated as one of the earliest diverging 
branches of the squamate tree, as an 
outgroup to iguanians in morphological 
phylogenies (Conrad, 2008; Gauthier et al., 
2012; but see Simdes et al., 2018). This 
results from Huehuecuetzpalli sharing some 
morphological characters with iguanians, 
such as a dorsal process of the squamosal 
and a reduced postfrontal located anteriorly 
to the postorbital (Fig. 2, chars. 1, 4). 
However, with the incorporation of molec- 
ular data, the combined evidence trees 
place iguanians within Toxicofera, thus 
“dragging” Huehuecuetzpalli with them 
inside Toxicofera and as the sister taxon to 
iguanians owing to their shared morpholog- 
ical attributes (Wiens et al., 2010; Pyron, 
2017: but see Reeder et al., 2015). In a 
similar manner, the morphological signal 
linking Mosasauria to snakes and angui- 
morphs drives mosasaurs to be nested 
within Toxicofera either as the sister taxon 
to snakes (Lee, 2009; Reeder et al., 2015; 
Pyron, 2017; Simées et al., 2018) or 
anguimorphs (Wiens et al., 2010; Simdes 
et al., 2020b). Finally, extinct subclades of 
crown families (e.g., glyptosaurines, a sub- 
clade of anguids) are found to be nested 
within anguids (inside Anguimorpha) and, 
therefore, inside Toxicofera (Wiens et al., 
2010; Reeder et al., 2015; Pyron, 2017). We 
predict that such extinct groups within 
major crown clades (e.g., glyptosaurines) 
are highly unlikely to have their phyloge- 
netic placement determined independently 


by the morphological signal in the presence 
of sufficiently large (e.g., ca. >8-10 kb) 
molecular sampling of extant species— 
except in the case of a rampant homopla- 
sy-driven morphological signal, such as the 
clustering of limb-reduced taxa. 

The major exception to the observed 
pattern described above is in the case of 
fossil taxa that do not have strong support 
from the morphological data to be the sister 
taxon to (or a subclade of) an extant clade. 
For instance, the controversial fossil lizard 
Sineoamphisbaena hexatabularis (Late Cre- 
taceous, Mongolia) has been previously 
proposed to be the sister lineage of amphis- 
baenians (see above; Wu et al., 1996: 
Gauthier et al., 2012), but other studies 
have placed it as a sister lineage to or a 
member of Borioteiioidea (Kearney, 2003; 
Lee, 2005b; Conrad, 2008) or estimated 
borioteiioids and Sineoamphisbaena as the 
successive sister lineages to amphisbaenians 
(Lee, 2009). The lack of consensus on the 
phylogenetic affinities of Sineoamphisbaena 
is further evidenced by the fact that all 
proposed systematic placements for this 
taxon are usually met with extremely low 
topological support (e.g., Bremer index of 1 
in maximum parsimony analyses). As such, 
the inclusion of molecular data and the 
nearly inevitable rearrangement of the 
relationships among extant families will 
likely affect and shift those controversial 
fossil taxa to potentially novel phylogenetic 
positions compared with the morphological 
tree (e.g., Lee, 2009; Wiens et al., 2010; 
Reeder et al., 2015; Pyron, 2017). However, 
as argued above, such rearrangement is 
most likely a result of their poor phyloge- 
netic signal rather than the actual capabil- 
ities of the morphological data to drive their 
phylogenetic placement independently of 
the molecular data. Such paradoxes can only 
be solved by improving character construc- 
tion and increasing character sampling in 
the situations where this is possible, given 
fragmentary remains. 
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Other fossil species with similarly non- 
consensual and poorly supported placement 
among different analyses include: Ardeo- 
saurus (Late Jurassic, Germany), the Early 
Cretaceous Calanguban, Tijubina, and 
Olindalacerta (Brazil); Tepexisaurus (Mex- 
ico); Scandensia, Pedrerasaurus, and Jucar- 
aseps (Spain); Dalinghosaurus (China); and 
Slavoia, Globaura, and Eoxanta (Late Cre- 
taceous, Mongolia) (Conrad, 2008; Bolet 
and Evans, 2012; Gauthier et al., 2012: 
Simoes et al., 2015a, 2017a). These species 
are also the ones that tend to reduce 
topological resolution and support (and 
most likely, tree accuracy) in combined 
evidence analyses and total evidence dating. 
We argue that algorithms for rogue taxon 
identification should be used to avoid wild- 
cards (see Pyron 2017), as evidence suggests 
they have a considerable negative effect for 
estimating accurate tree topologies under a 
variety of scenarios, both by increasing the 
overall amount of missing data and by 
destroying consensual support for optimal 
topologies (Luo et al., 2020; Vernygora et 
al., 2020). 


SYNTHESIS IN SQUAMATE EVOLUTION 
The Origin of Squamates 


Understanding the roots of squamate 
diversification, including the stem diver- 
gence from its sister clade (rhynchocepha- 
lians) and the placement of lepidosaurians 
in the greater context of diapsid reptile 
evolution, has remained elusive for nearly a 
century (Broom, 1925; Carroll, 1975, 1977). 
By far, the greatest limitations have been 
imposed by a patchy fossil record of 
squamates from the Triassic and Jurassic, 
mostly comprising fragmentary jaw ele- 
ments during the Jurassic (e.g., Waldman 
and Evans, 1994) anda complete absence of 
Triassic squamates until very recently. 
Another major factor was the lack of 
understanding of broadscale diapsid rela- 
tionships: Which taxa should act as the 
outgroup to lepidosaurs, and which clades 


make up Lepidosauromorpha? F inally, the 
overwhelming disagreement between most 
morphological and molecular hypotheses 
concerning the phylogeny of crown squa- 
mates (see above) has hindered a consen- 
sual approximation of the relationships 
between extant lineages, and thus our 
understanding of recent squamate evolu- 
tion. 

The recent discovery of the oldest known 
squamate from the Middle Triassic of Italy, 
followed by the first major agreement 
between morphological and molecular data- 
sets concerning important parts of the 
squamate tree, provides the first modern 
piece of evidence towards elucidating squa- 
mate origins (Simdes et al., 2018). By 
including a large number of other diapsid 
lineages, Simdes et al. (2018) determined 
that many diapsid reptile clades previously 
proposed as early-diverging lepidosauro- 
morphs, and thus supposedly the closest 
relatives to lepidosaurs (kuehneosaurids, 
younginiforms, sauropterygians, among oth- 
ers), actually fall in other parts of the 
diapsid Tree of Life. This latter conclusion 
was more thoroughly corroborated by the 
inclusion of all major clades of diapsids in a 
single data set by Simdes et al. (2018), 
although other studies with more limited 
taxon sampling already hinted towards this 
conclusion over the past 20 years. For 
instance, kuehneosaurids, gliding reptiles 
from the Triassic that were once thought to 
represent the earliest lizards (Robinson, 
1962; Estes, 1983), and later to be non- 
squamate lepidosauromorphs (Benton, 
1985; Gauthier et al., 1988a; Motani et al.. 
1998), have been recovered as non-lepido- 
sauromorphs in the previous large-scale 
phylogenies inclusive of most diapsid line- 
ages (Miiller, 2004; Hill, 2005: but see Chen 
et al., 2014). The same was detected for 
younginiforms, terrestrial and semiaquatic 
Paleozoic reptiles (Miiller, 2004: Hill, 2005), 
and more recently, for sauropterygians, 
semiaquatic and aquatic reptiles. (Hill, 
2005; Chen et al., 2014). 
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The growing current consensus is that the 
early part of the lepidosauromorph Tree of 
Life consisted of only a few stem taxa before 
the divergence between rhynchocephalians 
and squamates (e.g., Palaecagama, Sophine- 
ta, and Vellbergia, spanning the latest 
Permian [~252 Mal] to the early Middle 
Triassic [~230 Ma]) (Carroll, 1975, 1977: 
Evans and Borsuk-Bialynicka, 2009; Sobral 
et al., 2020). The identification of Mega- 
chirella (Figs. 1A-C) as the oldest stem 
squamate has now also been further cor- 
roborated by expansions and reanalyses of 
the largest available diapsid-squamate data 
set (Simoes et al., 2018: Garberoglio et al., 
2019a; Bittencourt et al., 2020: Simdes et 
al., 2020b: Sobral et al., 2020) by a new data 
set focused on sphenodontians (Simdes et 
al., 2020a) and by an independently built 
data set (Scheyer et al., 2020). Given that 
the oldest known fossil squamate (Mega- 
chirella) and rhynchocephalian (cf. Diphy- 
dontosaurus) come from the Middle 
Triassic at ~242 Ma and ~230 Ma, 
respectively, we are now able to narrow 
down the divergence between lizards and 
tuataras to the vicinity of the Permian— 
Triassic extinction. Indeed, this has been 
recovered by independent data sets imple- 
menting phylogenomic (Irisarri et al., 2017; 
Burbrink et al., 2020) and total evidence 
dating methods (Simées et al., 2018, 2020). 
Yet, it remains to be established whether 
the split between lizards and _ tuataras 
happened before or after the Permian— 
Triassic mass extinction. The latter will 
depend on even more precise estimates of 
divergence times and additional fossil ma- 
terial with the results being of great 
consequence for understanding the early 
drivers of squamate evolution. 

The radiation of crown squamates also 
finds consensual support on several grounds 
regarding the placement of gekkotans as the 
earliest evolving crown clade of squamates, 
along with dibamids (but see more below in 
Perspectives). Most notably, the early di- 
vergence of gekkotans in the squamate tree 
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has been recovered by essentially all mo- 
lecular data sets to date and finds agreement 
in the latest morphological hypotheses 
(Simoes et al., 2018; Simdes et al., 2020b). 
Importantly, this hypothesis also finds 
support in the fossil record, where the 
oldest known gekkotans come from the 
Middle Jurassic (Caldwell et al., 2015), 
and articulated stem gekkotans are known 
since the Late Jurassic (Simdes et al., 
2017a). In contrast, the oldest putative 
iguanian fossils are represented by fragmen- 
tary elements from the Early Cretaceous 
and more complete specimens only in the 
Late Cretaceous (Gao and Nessov, 1998: 
see also review by SimGes et al., 2017c). The 
fossil record may also support an earlier 
divergence of scincoids (the next major 
branch in the squamate tree following the 
consensual trees), because paramacellodids 
(Figs. 1F-I) also have a fossil record since 
the Middle Jurassic (Waldman and Evans, 
1994) and achieved a near cosmopolitan 
distribution by the Early Cretaceous (Bit- 
tencourt et al., 2020). 

Furthermore, it is important to recognize 
that the placement of gekkotans as part of 
the earliest radiation of crown squamates 
has substantial support from a morpholog- 
ical perspective. Gekkotans possess a vari- 
able number of features that are common to 
early evolving diapsid reptiles from the 
Triassic, including early lepidosaurs and 
most fossil sphenodontians, but which are 
absent in all other squamates. Those include 
the presence of amphicoelic vertebrae and 
the persistence of a notochordal canal in 
adults (both in gekkonines and diplodacty- 
lines, and also in stem squamates and 
Eichstaettisaurus), a perforated stapes (in 
eublepharines and in dibamids), the metotic 
fenestra undivided externally, although al- 
ready divided medially (gekkotans and 
xantusiids), and the presence of paired 
premaxillae (also in stem squamates and 
Eichstaettisaurus) (Simoes et al., 2018) (e.¢., 
Fig. 2, chars. 7, 8). Those anatomical 
attributes have historically been suggested 
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as the result of evolutionary reversals within 
gekkotans (Estes et al., 1988). However, as 
molecular hypotheses had already suggested 
(see above and Fig. 5), with further support 
by recent morphological phylogenies 
(Simoes et al., 2018; Simdes et al., 2020b), 
those attributes actually represent plesio- 
morphic traits in squamates. On the other 
hand, purported characters shared by igua- 
nians and sphenodontians are, in fact, the 
result of convergent evolution, which is 
further clarified by morphological aspects of 
stem squamates, especially Megachirella 
and Huehuecuetzpalli, that are shared with 
squamate outgroups and at least some 
gekkotans but are absent in all other 
squamates, such as the presence of amphi- 
coelic vertebrae and paired premaxillae. 
Therefore, the inclusion of closely related 
outgroup taxa to the early divergence of 
squamates in morphological phylogenies, 
such as early-evolving sphenodontians and 
early lepidosaurs, is fundamental to retrieve 
a more appropriate estimate of character 
evolution than the historically more com- 
mon and simplistic approach of including 
only the extant Sphenodon punctatus, which 
is separated from squamates by a branch of 
~250 million years. 


Phylogenetic Dating of Squamata 


Attempts to date the origin of major 
squamate lineages accompanied early esti- 
mates of the “molecular” hypothesis of 
squamate phylogeny (e.g., Vidal and Hedg- 
es, 2005). Those authors estimated a partic- 
ularly old age (240 Ma; 221-251 Ma range) 
for the MRCA of crown Squamata in the 
Middle Triassic. Subsequent authors began 
to refine these estimates (Mulcahy et al., 
2012; Jones et al., 2013) by integrating more 
complex models for molecular clock esti- 
mation and the incorporation of fossil data. 
These studies, and numerous subsequent 
estimates, have been remarkably consistent 
in estimating the MRCA of crown Squamata 
in the early Jurassic: ~180-190 Ma. The 


Timescale of Life database (TimeTree, 
2020) records at least 22 studies estimating 
the age of this node, ranging from 173 Ma 
(Gamble et al., 2015) to the 240 Ma of Vidal 
and Hedges (2005), with a median of 190 
Ma. However, these estimates are not 
independent of each other in terms of 
molecular or fossil data, relying as they do 
on similar loci and interpretations of existing 
fossil for use as node age constraints. Thus, 
these estimates are subject to the myriad 
concerns and considerations that affect 
secondary node age calibrations (Parham 
et al., 2012). 

These issues can be partially alleviated by 
total evidence dating that directly employs 
fossils as terminal taxa, incorporating mor- 
phological as well as molecular substitutions 
with direct observations of speciation, ex- 
tinction, and divergence times (Pyron, 2011; 
Ronquist et al., 2012). However, such 
methods still suffer from incomplete mod- 
eling of relevant processes and adequate 
associated priors (e.g., Arcila et al., 2015) 
and are thus still in their infancy. Only two 
sets of studies have attempted total evi- 
dence dating of Squamata to date, those of 
Pyron (2017) and Simoes et al. (2018) [and 
subsequent expansions of the latter (Gar- 
beroglio et al., 2019a; SimGes et al., 2020b)], 
incorporating both molecular and morpho- 
logical matrices (Table 1). 

Using the “standard” Markov K (Mk) 
models available for morphological evolu- 
tion under the fossilized birth—death pro- 
cess in MrBayes (see Ronquist et al., 2012), 
Pyron (2017) recovered an estimate of 190 
Ma for the MRCA of crown Squamata. This 
remarkable congruence with the DNA-only 
studies mentioned above is notable because 
he included only a single calibration on the 
root of Lepidosauria (238-250 Ma) from 
Jones et al. (2013). All other temporal 
calibration came from the observed ages of 
the fossil tips, with no internal node 
calibrations, and yet still recovered the 
same estimate of ~190 Ma as previous 
studies. Incorporating a more complex 
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Total Evidence Dating (Ma) 
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Clade Irisarri et al. (2017)4 


Burbrink et al. (2020)1 


Pyron (2017)} Simoes et al. (2020b)} 


Lepidosauria (crown ) 270.8 (246.2-298. 1) 241.5 244.9 (239.0-249.5) 259.6 (247.6-2,721.2) 
Squamata (total) = = = 250.3 (247.5-261.1) 
Squamata (crown) 203.7 (183.1—228) 193.2 186.6 (172.3-202.7) 186.5 (176.9-198.3) 


Gekkota (total) = 


159 (139.7-178.2) 154.2 (150-171.2) 


Gekkota (crown) 162.2 (124.7-196.2) 119.8 89 (72—103.9) 61.5 (50:2-75:4) 
Scincoidea (total) ~ 162.6 (153.1-173.6) 155.8 (150-168.9) 
Scincoidea (crown) 141.5 (99.9-175.2) 178.1 146.8 (131.6-161.5) 103.8 (79.5-127.8) 
Lacertoidea (total) =, = = = 
Lacertoidea (crown) 172.8 (112.4—200) ipzel 151.7 (136.1-164.9) 151 (139.5-164.8) 
Serpentes (total) = 130.3 (118.4-142.4) 166.6 (166.1-168)° 
Serpentes (crown) 69.5 (42.7-106.4) 129.8 121.9 (110.9-131.7) §3.6 (64.2—102.9) 
Anguimorpha (total) - ~ - 
Anguimorpha (crown) — 175.2 (155.8—198.1) 136.9 113 (101.4-126.1) 94.1 ($2.3-108.5) 


Iguania (total) 


Iguania (crown) 160.1 (142.9-183.6) 


157.8 


140.2 (127.9-151.5) 


129.2 (115.3-141.4) 102.5 (88.2121) 


+ Mean divergence times. Median values not made available in the original study. 

} Median divergence times (preferred relative to mean because it is less biased by outliers). 

° Minimum node age calibration applied to this clade on the basis of recent findings of the oldest snake remains from 
the Middle Jurassic (Caldwell et al., 2015). Pyron (2017) applied a tip-only dating calibration to this node, and so results 
are not directly comparable between studies. We also note that all phylogenomic studies are necessarily node-based only. 


model for morphological evolution (MkA), 
he estimated a slightly younger date of 187 
Ma, which is consistent with previous DNA- 
only studies showing that underparameter- 
ized models or undersampled data sets can 
inflate divergence time estimates (Schenk 
and Hufford, 2010; Mulcahy et al., 2012). 
Simées et al. (2018) employed essentially 
the same analytical strategy, but with a data 
set containing significantly more fossil taxa, 
particularly stem lineages. They tested 
calibrating the tree based on both tip-dated 
fossils only, as well as internal node plus tip 
age calibrations. They recovered a slightly 
older, Late Triassic age for crown Squamata 
at ~210 Ma, and ~257 Ma for the MRCA of 
all squamates (total group Squamata). Gar- 
beroglio et al. (2019a) essentially repeated 
the analyses of Simoes et al. (2018) with the 
addition of more fossil and extant snakes, 
yielding nearly identical results. However, 
in the latest expansion of this data set, 
Simoes et al. (2020b) employed approaches 


to avoid overestimating divergence times 
owing to deep root attraction in total 
evidence dating (Ronquist et al., 2016) with 
the use of two distinct software packages. 
The best performing model from SimG6es et 
al. (2020b) estimates times for the age of the 
MRCA of crown squamates congruent with 
the preferred analysis of Pyron (2017) 
implementing the Mka model at ~186.5 
Ma. Additionally, the age for the MRCA of 
all squamates was brought down to ~250 
Ma (immediately after the Permian—Triassic 
mass extinction) and the timing for the split 
between squamates and rhynchocephalians 
at ~259 Ma. 

Thus, the prospects for total evidence 
dating to converge widely accepted _esti- 
mates for squamate origins seem promising, 
especially when appropriate models are 
used to handle both morphological and 
molecular data. Chief among these are 
large-scale morphological matrices sampling 
more fossil taxa with well-constructed char- 
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acters that are comparable in extant species 
(e.g., Simées et al., 2017b). Another crucial 
avenue of investigation will be more accu- 
rately parameterized models for morpho- 
logical evolution (e.g., Harrison and 
Larsson, 2015; Bapst et al., 2016; Matzke 
and Wright, 2016; Wright et al., 2016; 
Pyron, 2017; SimGes et al., 2020a). Togeth- 
er, these total evidence dating approaches 
hold great promise for simultaneously 
resolving the topology of early squamate 
groups, as well as estimating their timescale. 


Toxicofera 


One of the most biologically intriguing 
developments in our recent understanding 
of squamate evolution has been the func- 
tional assemblage of salivary proteins (i.e., 
venoms) in Toxicofera. Human knowledge 
and fear of snake venoms is prehistoric and 
may have contributed to the evolution of 
primate vision systems (Isbell, 2006). 
Knowledge of venom in helodermatid 
lizards is much more recent (Woodson, 
1947). Anecdotally, belief that various other 
lizard species are venomous, including 
anguids, skinks, geckos, chameleons, and 
pleurodontans, is widespread in folklore 
around the world. Squamate venoms gen- 
erally consist of ordinary and abundant 
somatic proteins that undergo subsequent 
biochemical modifications that increase 
their toxicity (Hargreaves et al., 2014a). 
This process likely results either from 
genomic duplication and physiological re- 
cruitment from body tissues into modified 
salivary glands in the jaws or previous 
occurrence of these proteins in the salivary 
glands and subsequent restriction and 
localization into these glands therein. But 
how many times has this occurred? 

Recruitment of toxins was previously 
thought to represent a more recent evolu- 
tionary event restricted to the classically 
“venomous” snakes in the families Elapidae 
and Viperidae, with additional convergent 
origins in a few medically significant colu- 


brids. However, a series of biochemical 
analyses concurrent with the development 
of molecular hypotheses for squamate 
origins found widespread presence of “ven- 
om” proteins across snakes (Fry et al., 
2003). Subsequently, and contemporaneous 
with the expanded results of Vidal and 
Hedges (2005), Fry et al. (2006) announced 
the discovery of apparently functional ven- 
om genes, venom glands, and venoms across 
Anguimorpha, Iguania, and Serpentes, com- 
prising Toxicofera (gk. “those who bear 
toxins”), or the “venom clade.” 

This discovery implies a single early 
origin of “venom” in the MRCA of Tox- 
icofera ~170 Ma, with varying subsequent 
selection for toxicity and delivery apparatus. 
Although most of the >4,600 species in this 
clade are not medically significant to 
humans, bioactive salivary proteins may 
thus play a much larger role in prey capture 
than previously thought. For instance, the 
notably noxious and medically significant 
bite of the Komodo dragon (Varanus 
komodoensis) was long thought to be related 
to toxic oral bacteria obtained from eating 
carrion. It is now better understood to result 
(at least in part) from toxicoferan venoms 
(Fry et al., 2009), representing a second 
major lineage of medically significant lizards 
(Koludarov et al., 2017): even the acrodont 
iguanian Pogona, a common species in the 
pet trade, has a venom-secreting glandular 
apparatus in its jaws (Fry et al., 2006). In 
contrast, most other lizard species appear to 
have lost these apparatus and_ toxicity 
entirely. 

However, this interpretation of an early 
single origin of “venom” with a unique 
purpose for prey capture and a wide variety 
of subsequent evolutionary outcomes in 
terms of toxicity and complexity of appara- 
tus is not without controversy. Broader 
sampling of species and tissue types reveals 
widespread presence of many “venom”- 
type genes in nonfunctional roles across 
nontoxic structures (Hargreaves etival:, 
2014b). Even within snakes, many proposed 
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toxins are present across a variety of tissues. 
including salivary glands, and may thus have 
been restricted to those glands during 
independent origins of venom (Hargreaves 
et al., 2014a). There is also significant 
disagreement about the structure. function, 
and terminology of venoms and_ their 
apparatus within and among species (see 
Weinstein et al., 2012 and subsequent 
replies). 

Thus, it may be more accurate to think of 
the necessary genomic and_ physiological 
substrate for venoms (a particular suite of 
physiologically active proteins) to have 
evolved early in the history of Toxicofera. 
In contrast, the several lineages with 
generally agreed “venoms” and associated 
apparatus may represent multiple indepen- 
dent parallel origins from this shared 
substrate. Resolving this question is com- 
plex, requiring multiple lines of evidence, 
including genomics, biochemistry, physiolo- 
gy, histology, functional anatomy, phyloge- 
ny, and paleontological evidence. Well- 
defined ontologies regarding exactly what 
we consider “venom” and “venomous” are 
also a necessary complement to these data. 
Regardless, venom systems in several 
groups of toxicoferan squamates are among 
the most consequential and charismatic 
evolutionary innovations in the Tree of Life. 
A necessary component to this understand- 
ing, however, is an accurate resolution of 
the Squamate Tree of Life itself, with which 
to estimate the evolutionary history of these 
processes. 


The Origin of Snakes 


The placement of snakes within the 
squamate Tree of Life and the origins of 
the extreme morphological adaptations that 
characterize the snake body plan have long 
been a hot topic of debate among evolu- 
tionary biologists. Morphological data have 
long been considered limited in their 
capacity to understand the phylogenetic 
origins of snakes, owing to the rampant 


occurrence of apparent homoplasies uniting 
limb-reduced and skull-miniaturized taxa. 
discussed in detail above. Although molec- 
ular data have provided considerable sup- 
port for the placement of snakes with 
iguanians and anguimorphans within Tox- 
icofera, they have not been able to indicate 
their exact systematic placement within the 
latter. 

Despite limitations of their own, many 
studies combining morphological and mo- 
lecular data with extensive sampling of 
extant and fossil taxa have placed snakes 
within Toxicofera, and as the sister clade to 
Mosasauria, forming the clade Pythonomor- 
pha (Reeder et al., 2015; Pyron, 2017; 
Sim6des et al., 2018; Garberoglio et al., 
2019a). Accordingly, Pythonomorpha is the 
sister clade to Iguania+Anguimorpha 
(Reeder et al., 2015; Simdes et al., 2018: 
Garberoglio et al., 2019a). Given the 
ambiguous molecular support for the place- 
ment of snakes within toxicoferans, mor- 
phological data thus play a considerable role 
in inducing the latter configuration. As such, 
the choice of the morphological matrix will 
affect the placement of snakes within 
toxicoferans, and variations on the latter 
pattern may occur (e.g., Pyron, 2017). 
Therefore, the continual assessment and 
reevaluation of morphological characters 
and fossil specimens remains critical. 

In the last decade, important advances 
have been made for understanding the 
homology of components of the snake 
phenotype, aspects that have historically 
had an important effect in the placement of 
snakes on the basis of morphological data. 
Those aspects include: homology of the 
“postorbital element,” now recognized as 
the jugal bone instead of a postorbital or 
postfrontal (Palci and Caldwell, 2013; Gar- 
beroglio et al., 2019a); homology of the 
cervical and caudal peduncles, now recog- 
nized as homologous to the hypapophyses 
and haemapophyses of other squamates; 
and the absence of a crista circumfenestralis 
enclosing the stapes in the earliest snake 
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fossils (Palei and Caldwell, 2013; Garber- 
oglio et al., 2019a,b). New and more 
complete fossils of Cretaceous snakes have 
also continually expanded our understand- 
ing of the acquisition of various features 
characteristic of the body plan of crown 
snakes (Garberoglio et al., 2019a). Those 
findings have provided increasing morpho- 
logical support for the (Pythonomorpha 
[Iguania+Anguimorpha]) configuration of 
Toxicofera. 

However, the sister-group relationship 
between snakes and mosasaurs estimated 
by many recent analyses creates an addi- 
tional problem in squamate evolution: a 
ghost lineage of nearly 40 million years 
exists between the snake—mosasaur split at 
~171 Ma (SimGes et al., 2018; Garberoglio 
et al., 2019a) and the earliest mosasaurians 
in the fossil record ~132.9-125 Ma (Camp- 
bell Mekarski et al., 2019). This gap suggests 
that either further rearrangements of this 
topological configuration will emerge on the 
basis of new fossil and morphological data, 
or that a major gap in the early fossil record 
of mosasaurians remains to be filled. Simdes 
et al. (2020b) recently applied diversifica- 
tion priors that penalize excessively long 
ghost lineages (which negatively affect 
divergence time estimates), recovering mo- 
sasaurians as the sister lineage to angui- 
morphs instead of snakes, despite a strong 
morphological signal linking mosasaurians 
and snakes. This investigation indicated that 
available models are capable of efficiently 
penalizing such long ghost lineages, despite 
a strong signal in the data (see also Pyron, 
2017). However, because the earliest mem- 
bers of the mosasaurian lineage were, at 
some point, most likely inhabiting a terres- 
trial environment, it is possible that such 
long ghost lineage could be the result of the 
poor fossil record of terrestrial squamates 
during most of the Mesozoic until the late 
Karly Cretaceous (Cleary et al., 2018; Close 
et al., 2019). Further sampling of Late 
Jurassic and Early Cretaceous squamates is 
essential to elucidate this matter. 


Despite an apparent trend towards the 
consensual placement of snakes in the 
Squamate Tree of Life, those new results 
provide little insight regarding the environ- 
mental drivers of the early evolution of the 
snake body plan. The long and heated 
debate between a fossorial versus aquatic 
hypotheses for the origin of snakes (Lee, 
1998; Rieppel and Zaher, 2000; Lee, 
2005a,b; Gauthier et al., 2012; Hsiang et 
al., 2015) inherently depends not only on 
understanding the sister group relationships 
between early snakes and other squamates, 
but also a good understanding of the 
environments where the earliest snakes first 
evolved and their community structure. 
This knowledge is currently lacking, in that 
the oldest known snake fragments come 
from a variety of environments such as coal 
swamps, mixed coastal lake and pond 
systems, and epicontinental islands (Cald- 
well et al., 2015), and little is known of their 
ecological context and co-occurrence with 
other ecologically relevant species. If con- 
sensual answers are ever to emerge to settle 
this debate, they will depend on an 
increased knowledge of the fossil record 
and paleoenvironments of mid-Jurassic 
snakes and other squamates, sources of data 
outside the strict confines of molecular or 
morphological systematics. 

The phylogenetic origin of snakes has also 
stirred broad avenues of research towards 
understanding the origins of the extreme 
adaptations characterizing the snake body 
plan. Although limb reduction is widespread 
among squamates (Wiens et al., 2006), most 
groups retain at least partial development of 
their limb skeletons, with the major excep- 
tion being snakes. Although some relatively 
early diverging snakes retain remnants of 
their pectoral girdle (e.g., “henophidians”, 
such as boas and pythons), the vast majority 
of snakes have lost their limbs entirely. — 

Recently, genomic studies have revealed 
fundamental clues about the genetic mech- 
anisms responsible for the loss of limbs in 
snakes. The loss of limbs in Squamates 
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leading to the origin of snakes has been 
demonstrated to stem from mutations on an 
enhancer—the zone of polarizing activity 
regulatory sequence (ZRS)—of the sonic 
hedgehog gene (Shh). It was already known 
that this was a conserved regulatory region 
(Sagai et al., 2004) and that mutations on 
this enhancer can result into limb malfor- 
mation in vertebrates (Sagai et al., 2005; Hill 
and Lettice, 2013), and so it would be 
expected that snakes would have a high 
number of mutations on the ZRS. Surpris- 
ingly, however, this enhancer has been 
observed to be extremely conserved across 
extant early-evolving snakes retaining rem- 
nants of the hindlimbs (e.g., boas and 
pythons) (Kvon et al., 2016). Only later- 
evolving snakes without limb remnants 
demonstrate a much higher degree of 
substitutions on this enhancer, suggesting 
a progressive loss of function of the ZRS in 
snake evolution (Kvon et al., 2016). This 
hypothesis of a more gradual evolutionary 
loss of limbs in snakes has been later 
supported by new fossil evidence and 
phylogenetic reconstructions suggesting 
that early snakes retained relatively well- 
developed hindlimbs for at least the first 80 
million years since their origin (Garberoglio 
et al., 2019a). Combined, these results 
suggests that early snake evolution is not 
characterized by drastic changes on their 
limbs, but rather, on the skull (see also 
Caldwell, 2019). This idea finds further 
support from recently detected high rates 
of morphological evolution of the snake 
skull among the earliest branches of the 
snake evolutionary tree (Watanabe et al., 
2019) and early snake fossils already pos- 
sessing mandibular features extremely very 
similar to modern day taxa (Caldwell et al., 


2015). 


FUTURE PERSPECTIVES 
Phylogenetic Signal and Fossil Data 


Despite the long-term debate between 
the role of morphology and molecules in 


systematics (e.g., Scotland et al., 2003; 
Wiens, 2004), this issue was mostly con- 
cerned with the quality of such data sets per 
se for systematic purposes (regardless of 
analytical procedure) or focused only on 
more traditional analytical approaches such 
as maximum parsimony. However, the 
increasing availability of phylogenomic data 
sets and the emergence of total evidence 
dating with relaxed clock Bayesian inference 
under rapidly improving analytical proce- 
dures, makes it important to assess the 
effect of and the best approach towards 
combining morphological (from fossils and 
extant taxa) and molecular data to recon- 
struct the squamate Tree of Life. 

As illustrated in the sections above, 
evidence is ample for the limited role of 
morphological data and fossil taxa in driving 
estimated topological relationships in squa- 
mates when combined with moderate to 
large-sized molecular data—for either non- 
clock or relaxed clock total evidence dating. 
This situation may only be significantly 
altered by character weighting to balance 
out the overwhelmingly large number of 
molecular characters (Giribet, 2010), which 
has thus far not been attempted in large 
squamate data sets. However, even in the 
absence of character weighting, we argue 
that, despite such limitations, morphological 
data still play an important role for the 
placement of both extant and extinct taxa in 
the instances in which the molecular signal 
is weak or ambiguous (see Wiens et al., 
2010; Pyron, 2017). A clear example comes 
from the ambiguous molecular signal de- 
scribed above concerning the earliest di- 
chotomy in the squamate tree (dibamids 
and gekkotans) and over the relationship of 
the three major toxicoferan clades. In such 
situations, the morphological signal drives 
the placement of gekkotans as the earliest 
evolving crown squamate lineage relative to 
dibamids on the basis of the plesiomorphic 
features of some gekkotans in the context of 
early diapsid and early lepidosaur evolution, 
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and which are not present (or inapplicable) 
in dibamids (see above; Simoes et al., 2018). 

This signal is so strong that even other 
morphological data sets without a broad 
sampling of early diapsids and lepidosaurs 
(including only a handful of sphenodontian 
taxa for outgroup comparison) also recover 
gekkotans as diverging earlier relative to 
dibamids when combined with molecular 
data, using both the Conrad (2008) and 
Gauthier et al. (2012) data sets tested by 
Pyron (2017). The morphological signal may 
also play an important role in the placement 
of snakes within Toxicofera given the 
ambiguous molecular signal, although this 
has yet to be formally tested. F urthermore, 
morphological data play an important role 
even in the presence of extremely large 
molecular data sets by estimating the 
relationships within entirely extinct lineag- 
es, such as mosasaurians, borioteiioids, and 
priscagamids, among many others, which 
will also affect their sister group relation- 
ships to extant clades (see Wiens et al., 
2010; Pyron, 2017). Finally, the inclusion of 
both morphological and molecular data has 
a fundamental role in total evidence dating, 
by avoiding ad hoc hypotheses for the 
placement of fossil taxa and further avoiding 
a priori constraining of topological relation- 
ships (Ronquist et al., 2012; Zhang et al., 
2016). 

Another important point of consideration 
is that a great portion of the Tree of Life, 
especially in early-diverging branches, is 
represented only by fossils, with some 
accounts estimating that up to 99.9% of all 
species that ever lived are now extinct 
(Novacek and Wheeler, 1992). This obser- 
vation implies an inherent property of 
phylogenetic trees of clades with extant 
representatives: The broader the taxonomic 
sampling, or the deeper in time we need to 
go to find the last common ancestor of all 
sampled species, the greater the likelihood 
of encountering entirely extinct lineages. As 
noted above, the phylogenetic relationships 
of these extinct lineages can thus only be 


informed by morphological (usually osteo- 
logical) characters and data sets. 

A fundamental corollary of this property 
is the prediction that the broader or deeper 
in time the taxonomic sampling (e.g., higher 
level phylogenies) extends, the greater the 
chances of creating artifactual sampling 
gaps if fossils are not included (e.g., 
compare Figs. 5 and 6). We call this 
property the “deep-time sampling bias.” 
Thus, many, if not most, “long branches” in 
extant-only trees (e.g., the branch leading to 
Sphenodon) are not actually long; they are 
simply undersampled, creating increasingly 
larger taxon sampling gaps irrespective of 
the sampling effort of extant taxa. Concom- 
itantly, “breaking” these long branches by 
sampling additional fossil taxa is well known 
to exert a strong and presumably positive 
effect on phylogenetic inference and dating 
(Gauthier et al., 1988b; Wiens, 2004; 
Guillerme and Cooper, 2016; Luo et al., 
2020). 

In the case of squamates, the lineages 
representing most extant genera and fami- 
lies coalesce at some point in the last 75- 
100 million years (Pyron, 2017; SimGes et 
al., 2018; Burbrink et al., 2020; Simdes et 
al., 2020b), but the divergence between the 
major squamate higher clades is much more 
ancient, concentrated between 75 and 175 
Ma (Table 1; Pyron, 2017; SimGes et al., 
2018; Burbrink et al., 2020; Simdes et al., 
2020b), sometimes marked by some long 
stem branches. As expected, the 100 million 
years between the origin of the higher 
squamate clades and the divergence of most 
extant families represent a timeframe 
densely interpermeated by extinct families 
(e.g., Mosasauria, Borioteiioidea, Parama- 
cellodidae; Figs. 1J—N) and also many stem 
members of extant lineages (e.g., Eichstaet- 
tisaurus and Gobekko for gekkotans; Najash 
and Dinilysia for snakes; Figs. 1D, E). 
Fossil squamates ranging in age from ~75 
Ma to the time of the oldest known fossil 
squamate (242 Ma) are thus of particular 
relevance to sample several extinct lineages 
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and break ancient long branches and 
avoiding deep fossil sampling biases. 

Therefore, genomic data are most likely 
to provide highly accurate and robust 
phylogenetic hypothesis relationships 
among younger or more recent lineages 
given the recent radiation of the sampled 
species (e.g., family or subfamily level). 
However, higher level phylogenies have a 
much greater need for sampling of extinct 
lineages for a balanced taxonomic represen- 
tativeness, thus necessitating a morpholog- 
ical partition for the inclusion of fossil taxa. 
Recent studies demonstrate the latter to be 
fundamental for improving accuracy for 
divergence time estimation (Guillerme and 
Cooper, 2016; Luo et al., 2020). Although 
such an effect was already known in some 
sense (e.g., Graybeal, 1998; Wiens, 2003), 
the particular example of the squamate Tree 
of Life throws it into sharp relief for 
broadscale phylogenetic estimation of tem- 
porally diverse, ancient groups. 


Whole Genomes 


Genomics is still in relative infancy for 
evolutionary and phylogenetic studies, and 
whole-genome data present a wealth of data 
for resolving relationships among lineages 
and clarifying the evolutionary origin of 
traits that have not yet been fully deployed 
for these purposes, particularly in squa- 
mates. For instance, analysis of gene order 
conservation (synteny mapping) represents 
a promising source of phylogenetic signal— 
essentially measuring the “morphology” of 
the genome. Analyses thereof have resolved 
contentious branches in prokaryotes (Shif- 
man et al., 2014) and mammals (Luo et al., 
2012). Other sophisticated methods of 
analyzing whole-genome data include topol- 
ogy-by-location scans, wherein genealogies 
are estimated for loci across the whole 
genome, showing complete patterns of gene 
tree variation (Li et al., 2019). 

One can thereby highlight the empirical 
distribution of discordant gene trees across 


areas of the genome that are duplicated, 
introgressed, or vary in effective population 
size, giving more decisive estimates of the 
overall species tree (e.g., Ravinet et al., 
2018). Such approaches might be particu- 
larly applicable to longer internal or termi- 
nal branches, such as the placement of 
Dibamia or amphisbaenian monophyly, by 
revealing which topologies are linked most 
closely to highly conserved regions of the 
genome (see Singhal et al., 2021). 

Other sources of phylogenetic signal 
abundant in whole-genome data are mobile 
genetic elements such as SINEs and 
LINEs, which have long been noted for 
their potential applicability to difficult 
phylogenetic problems (Hillis, 1999), but 
rarely deployed at scale. Genomic scans can 
produce a wealth of such data as bycatch 
from other types of loci, offering a corrob- 
orative counterpoint to gene tree heteroge- 
neity because of short internal branches and 
incomplete lineage sorting (Cloutier et al., 
2019). Such approaches may be confounded 
by the complex diversity and history of such 
elements in squamate genomes (Pasquesi et 
al., 2018). Nevertheless, element-based 
approaches might be helpful in discriminat- 
ing among topologies in the anomaly zone of 
rapid radiations, such as Pleurodonta and 
booid, colubroid, and elapoid snakes. 

Genomic data are, of course, also indis- 
pensable in understanding the origin, evo- 
lution, and development of complex traits, 
in conjunction with detailed morphological 
and physiological observations. Myriad ex- 
amples are present in the recent literature, 
such as the origins of visual and olfactory 
adaptations for prey capture (Perry et al., 
2018), venom evolution (Schield et al., 
2019), metabolic adaptations (Lind et al., 
2019), sex chromosomes (Gamble et. al. 
2017), and viviparity (Gao et al. 2019). 
Genomic assemblies of squamate species 
are becoming increasingly common. for 
snakes (Kerkkamp et al., 2016) and anno- 
tated drafts available for some lizards 
species within iguanians (Alfoldi et. al., 
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2011), anguimorphans (Song et ale, QOMS), 
and geckos (Xiong et al., 2016). As of 28 
May 2020, more than 241 genome-scale 
data sets are present on GenBank, providing 
a wealth of data for both phylogenetic 
inference and comparative analysis. 


CONCLUSIONS 


Understanding the overall relationships 
among the numerous lineages of living and 
extinct squamates has been one of the major 
goals in herpetology and evolutionary biol- 
ogy for over a century. This long-lasting 
problem has generated fierce debates across 
disciplines towards understanding some of 
the most fundamental questions regarding 
the evolution of lizards, snakes, amphisbae- 
nians, and tuataras. However, the last two 
decades have been especially promising 
toward clarifying some of the most disputed 
aspects of squamate evolution. 

Squamate phylogenomics has revealed 
that nearly all genomic regions result in 
very similar species trees, with only a few 
points of standing disagreement. Important- 
ly, agreement between morphological and 
molecular hypotheses is starting to emerge. 
Those hypotheses, along with combining 
morphological and molecular data in total 
evidence dating, have enabled the develop- 
ment of important points of increasingly 
consensual agreement and synthesis in 
squamate evolution. We propose that major 
points of synthesis in squamate evolution 
are represented by a much improved 
understanding of the early fossil record 
and origin of squamates from other reptile 
lineages—a remarkable convergence of 
estimates of diversification times for the 
major lineages of squamates from. total 
evidence dating, monophyly of Toxicofera, 
and the placement of snakes as the sister 
lineage to iguanians+anguimorphans, irre- 
spective of the placement of mosasaurians. 
Other emerging patterns, but which need 
further corroboration, include monophyly of 


amphisbaenians and paraphyly of scoleco- 
phidians. 

The major standing problems in squa- 
mate phylogenetics that remain to be solved 
by larger scale phylogenomics and morpho- 
logical data sets (most crucially, including 
more fossil lineages) are the internal rela- 
tionships of pleurodont iguanians and co- 
lubroid snakes, the placement of dibamids, 
and the placement of many fossil families 
(e.g., mosasaurians and paramacellodids). 
The relationships of iguanians, colubroids, 
and dibamids will be resolved, one hopes, 
by genome-scale data sets, as well as 
possibly by the discovery and inclusion of 
new fossil taxa in morphological matrices. 
Estimating fossil relationships and more 
accurately dating and reconstructing deep 
branches in the squamate tree depend on 
expanding our currently scarce knowledge 
of the squamate fossil record during the first 
half of squamate evolution: from the Middle 
Triassic to the Early Cretaceous. Perhaps 
more than any other factor, this is the most 
essential source of data to understand not 
only early squamate taxonomy, but also to 
improve our ability to break long stem 
branches and reconstruct early patterns of 
squamate evolution and the time of origin of 
its major clades. 

As evidence accumulates regarding the 
benefits of combining morphological and 
molecular data for inferring accurate time- 
trees, the future of estimating higher level 
evolutionary trees seems to be inherently 
dependent on accurate estimates of both 
morphological and molecular evolution. 
Therefore, efforts should be put into 
improving models and protocols to assem- 
ble and analyze such data sets, both 
independently as well as combined in total 
evidence dating. Finding agreement be- 
tween both data types is an important 
starting point. It may be impossible to 
obtain a full agreement between morpho- 
logical and molecular trees, because they 
represent fundamentally distinct dimen- 
sions of organismal evolution and are 
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differentially affected by selective pressures 
and by fixed and plastic variation in 
development and ontogeny. However, at 
least partial levels of agreement are achiev- 
able and necessary to avoid conflicting 
phylogenetic signals that negatively affect 
timetree inference. The development of 
more complex models of morphological 
character evolution is paramount to a more 
biologically realistic modeling of the phe- 
notype, comparable to the array of molec- 
ular substitution models currently available. 
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